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PREFACE 


The  elastic  behavior  of  composite  materials  is  a  subject  of 
strong  topical  interest  because  of  the  technological  applications 
of  these  materials.  Since  no  real  solid  is  infinite,  it  is 
important  to  study  the  effect  of  a  free  surface  on  the  elastic 
properties  of  composite  materials.  A  topic  of  special  interest  is 
the  stress  distribution  due  to  an  applied  load  in  such  solids. 
Many  papers  have  already  been  published  on  this  topic.  The 
behavior  of  stress  near  the  intersection  of  the  free  surface  and 
the  interface  in  a  bimaterial  composite  is  singular  and  quite 
complex  (see,  for  example,  [1-14]  and  other  papers  quoted  in  these 
references) .  This  is  usually  referred  to  as  the  free-edge  effect 
in  composite  materials. 

The  analysis  of  the  free-edge  effect  is  particularly 
important  because  damage  in  a  composite  material  initiates  at  the 
free-edges,  often  as  ply  cracks  and  delaminations .  Usually,  these 
delaminations  are  the  source  of  compressive  failure  in  reverse 
fatigue  loading  [15].  The  free-edge  effect  is  a  characteristic 
property  of  composite  materials  and,  in  general,  not  exhibited  by 
homogeneous  solids.  Such  a  distinction  between  homogeneous  and 
composite  solids  should  not  be  surprising  because  a  similar 
distinction  is  also  found  in  the  behavior  of  the  stress  near  an 
interfacial  crack  tip.  Unlike  that  in  a  homogeneous  solid,  the 
stress  near  a  crack  tip,   shows  strong  oscillations  [16,17]. 

There  is  some  confusion  in  the  literature  regarding  the 
nature  of  singularities  in  the  stress  field  near  the  free  surface. 
In  general  it  has  been  found  by  using  different  numerical 
techniques  that  the  stress  behaves  like  r  near  r=0  where  r  is  the 
radius  vector  and  8  is  a  negative  number — usually  a  fraction. 
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Zwiers  et  al.  [10]  have  shown  that  in  case  of  an  out-of-plane 
load  (the  generalized  plane-strain  problem) ,  logarithmic 
singularities  are  also  present  in  the  stress  field.  Further,  it 
has  been  pointed  out  [10]  that  additional  singularities  like 
[log(r)]n  (n-1, 2 , 3 , . . . )  may  also  arise  in  the  stress  field. 
Moreover,  the  behavior  of  the  nonsingular  terms  has  not  been 
properly  investigated  so  far  except  numerically.  There  is  clearly 
a  need  for  a  reliable  complete  solution  which  gives  the  stress 
distribution  in  the  composite  solid  containing  a  free  surface. 

The  stress  distribution  as  well  as  many  other  elastic 
properties  of  a  solid  can  be  easily  calculated  by  using  the 
elastic  Green's  function  of  the  solid.  The  Green's  function  gives 
the  solution  of  the  elastic  equilibrium  equations  for  a  unit  force 
subject  to  all  the  prescribed  boundary  conditions.  Once  the 
Green's  function  is  known,  the  solution  of  the  elastic  equilibrium 
equations  can  then  be  obtained  for  any  arbitrary  force 
distribution  by  a  straightforward  integration.  The  Green's 
function  is  therefore  a  very  powerful  mathematical  tool  for 
calculating  the  elastic  response  of  a  solid. 

The  objective  of  the  present  work  is  to  calculate  the  elastic 
Green's  function  for  an  anisotropic  composite  solid  containing  a 
free  surface  and  apply  it  to  calculate  the  stress  distribution  in 
the  solid  to  an  applied  load.  We  shall  consider  a  semi-infinite 
bimaterial  composite  solid  containing  a  planar  interface  and  a 
single  free  surface  normal  to  the  interface.  The  solid  is  assumed 
to  extend  to  infinity  in  other  directions.  This  problem  will  be 
referred  to  mathematically  as  the  surface-interface  problem  and 
the  corresponding  Green's  function  as  the  surface-interface 
(abbreviated  as  the  SI)  Green's  function. 

We  shall  calculate  the  SI  Green's  functions  for  the 
plane-strain  and  the  generalized  plane-strain  problems.  The 
out-of-plane    components    of    the    strain    tensor    are    zero    in  the 
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plane-strain  case  and  are  nonzero  but  constant  in  the  generalized 
plane-strain  case.  The  out-of -plane  components  of  the  displacement 
and  the  stress  tensor  will  be,   in  general,  nonzero,   in  both  cases. 

The  SI  Green's  function  as  calculated  here  is  applicable  to 
different  types  of  bimaterial  composites  such  as,  for  example, 
those  containing  different  materials  separated  by  a  phase  boundary 
or  different  orientations  of  the  same  crystal  separated  by  a  grain 
boundary.  It  can  also  be  applied  to  laminated  fiber  composites  in 
which  different  layers  have  different  fiber  orientations. 

Normally,  the  Green's  function  method  is  used  to  calculate 
the  stress  distribution  in  a  solid  due  to  an  specified  force.  We 
have,  however,  extended  the  Green's  function  method  to  the  case 
when  instead  of  a  force,  an  out-of-plane  displacement  field  is 
specified.  This  enables  us  to  study  the  technologically  relevant 
generalized  plane-strain  problems. 

In  two  earlier  papers  [18,19],  henceforth  referred  to  as  I 
and  II,  we  have  calculated  the  elastic  plane-strain  Green's 
function  for  an  infinite  bimaterial  composite  solid.  This  Green's 
function  gives  the  solution  of  the  elastic  equilibrium  equations 
and  satisfies  the  continuity  conditions  at  the  interface.  In  this 
paper  we  shall  use  this  Green's  function  to  calculate  the  SI 
Green's  function  which,  in  addition,  will  also  satisfy  the  free 
surface  condition  at  the  free  surface. 

The  calculation  of  the  SI  Green's  function  would  require  the 
solution  of  an  inhomogeneous  generalized  vector  Hilbert  problem. 
The  properties  and  the  solution  of  the  corresponding  scalar 
problem  have  been  discussed  in  the  excellent  treatise  by 
Muskhelishvili  [20]  and  the  vector  problem  by  Vekua  [21]  and  also 
briefly  in  [20]  .  However,  in  the  present  case,  we  find  it  more 
convenient  to  solve  the  aforementioned  Hilbert  problem  by  an 
alternative    method    using    a    complex    transform.     This    gives  the 
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result  in  a  closed  integral  form  which  can  be  calculated  quite 
easily.  This  method  should  be  generally  useful  in  several  problems 
in  stress  analysis  in  composite  materials  containing  cracks, 
interfaces,  and  free  surfaces. 

The  earlier  methods  [1-14]  for  calculation  of  the  stress 
distribution  near  a  free  surface  in  a  laminated  composite  are 
either  purely  numerical  or  use  a  series  representation.  The  series 
representation  is  valid  only  near  the  origin.  Zwiers  et  al.  [10] 
have  used  the  Stroh's  method  [22]  for  obtaining  a  series 
representation  of  the  stress  in  terms  of  the  roots  of  a  12  x  12 
complex  matrix.  All  these  calculations  were  done  for  a  solid 
subjected  to  certain  specific  loadings. 

The  purely  numerical  methods  are  not  reliable  near  the  origin 
where  the  stress  is  singular  and  quite  complex.  As  pointed  out  in 
[10],  a  series  representation  obtained  by  using  an  eigenf unction 
expansion  method  [7]  or  Lekhnitskii  method  suffers  from  lack  of 
completeness  of  eigenfunction  expansion  as  well  as  numerical 
convergence.  Moreover,  some  of  these  methods  do  not  account  for 
out-of-plane  components  of  the  displacement  and  the  stress  tensor 
which,   in  general,  will  exist  even  in  plane-strain  problems. 

The  Green's  function  method  can  be  applied  to  solids 
subjected  to  any  loading  or  force  distribution.  In  addition  to  the 
stress  distribution,  it  can  be  used  to  calculate  the  interaction 
between  various  defects  such  as  dislocations,  surfaces,  and 
cracks.  The  interaction  between  defects  determines  several 
important  physical  properties  of  materials  such  as,  for  example, 
the  work  hardening  of  materials. 

The  Green's  function  method  accounts  for  the  out-of-plane 
components  of  the  displacement  field  and  the  stress  tensor  in  the 
plane-strain  problems.    It  gives   the  result   in  a   closed  integral 
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representation  which  can  be  evaluated  numerically  quite  easily. 
The  result  is  therefore  valid  for  the  entire  solid  and  not  only 
near  the  origin.  The  integral  for  the  Green's  function  can  also  be 
evaluated  analytically,  which  is  useful  in  identifying  the  exact 
nature  of  all  the  singularities  which  may  be  present  in  the 
stress.  This  method  is  also  numerically  more  convenient  than  other 
methods  because  it  involves  manipulation  of  only  a  6  x  6  complex 
matrix  and  not  12  x  12 . 

This  paper  is  in  two  parts.  In  Part  I,  we  have  calculated  the 
SI  Green's  function  and  applied  it  to  plane-strain  problems  in 
cubic  solids.  Specifically,  we  have  calculated  the  stress 
distribution  due  to  a  line  load  in  a  cubic  solid  containing  a  Z  5 
grain  boundary  and  a  free  surface  normal  to  the  grain  boundary 
interface.  Both  the  plane-strain  and  antiplane-strain  problems 
have  been  considered.  In  Part  II,  we  have  described  an  extension 
of  the  Green's  function  method  to  the  generalized-plane-strain 
problem  which  enables  us  to  calculate  the  elastic  response  of  a 
solid  to  a  prescribed  out-of -plane  strain.  We  have  applied  it  to 
calculate  the  stress  distribution  in  a  composite  solid  subjected 
to  an  out-of -plane  load. 

Our  final  results  are  in  the  form  of  an  analytical  closed 
integral  representation  which  can  be  calculated  numerically  as 
well  as  analytically.  The  analytical  evaluation  of  the  integral  by 
using  the  contour  integration  method  gives  the  result  in  the  form 
of  a  series.  This  series  contains  singular  as  well  as  non-singular 
terms.  In  general,  it  is  difficult  to  ascertain  the  convergence  of 
the  series  representation.  The  numerical  evaluation  of  the 
integral  does  not  suffer  from  such  problems  and  gives  the  result 
to  arbitrary  accuracy. 

The  analytical  series  representation  can  be  used  to  identify 
precisely  the  nature  and  the  weight  of  each  singularity  in  the 
stress    distribution.    This    is    useful    in    developing    a  numerical 
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finite-element  technique  for  the  stress  analysis  in  a  real  solid. 
A  purely  numerical  finite-element  technique,  when  applied  to 
composite  materials  containing  free  surfaces,  suffers  from  poor 
convergence  near  the  origin  because  of  the  singularities  in  the 
stress.  A  knowledge  of  the  analytical  behavior  of  the  stress  near 
the  origin  can  be  used  to  generate  the  so-called  enriched  elements 
near  the  origin,  which  substantially  improves  the  convergence  of 
the  numerical  technique  [14]. 
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ABSTRACT 


This  paper  is  in  two  parts.  In  part  I,  the  elastic 
plane-strain  Green's  function  is  calculated  for  an  anisotropic 
bimaterial  composite  solid  containing  a  free  surface  normal  to  the 
interface.  An  exact  integral  representation  is  obtained  for  the 
Green's  function  which  is  evaluated  numerically.  The  integral  is 
also  evaluated  analytically,  which  gives  a  series  representation 
for  the  Green's  function.  For  illustration,  the  formalism  is 
applied  .  to  the  antiplane-strain  problem  and  the  plane-strain 
problem  in  a  cubic  solid  containing  a  Z-5  grain  boundary.  In  part 
II,  the  Green's  function  is  extended  to  the  case  of  generalized 
plane  strain  and  applied  to  calculate  the  stress  field  in  a  cubic 
solid  containing  a  Z  5  boundary  and  to  fiber-reinforced  composites 
with  various  fiber  orientations  subjected  to  an  out-of -plane  load. 

It  is  found  that,  as  predicted  by  earlier  authors,  the  stress 

is    singular    at    the    intersection    of    the    free    surface    and  the 

interface  in  both  plane-strain  and  generalized  plane-strain  case. 

These    singularities     in    the    stress    field    associated    with  the 

presence   of   the   free   surface   are   identified   and   discussed.  The 

singularities   in  the   stress   field  are   found  to   be   of      the  type 
—5 

r  ,  ln(r)  as  well  as  containing  higher  powers  of  ln(r),  where  5 
is  between  0  and  1  and  r  is  the  radial  distance  from  the 
intersection  of  the  free  surface  and  the  interface.  In  addition, 
it  is  found  that,  in  general,  the  stress  field  will  also  contain 
an  oscillatory  factor  of  the  type  exp [ igln (r) ] ,  where  g  is  a 
constant  which  depends  upon  the  material  parameters  of  the  two 
solids.  However,  for  the  case  of  the  generalized  plane  strain  in 
fiber-reinforced  composites,  the  factor  g  is  almost  0  so  that  the 
oscillatory  behavior  is  absent. 

Part  I  also  gives  a  simple  and  convenient  method  for  solving 
an  inhomogeneous  generalized  vector  Hilbert  problem  which  contains 
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a  nonsingular  kernel  in  addition  to  the  usual  singular  Hilbert 
kernel.  The  solution  of  this  Hilbert  problem  is  required  for 
obtaining  the  Green's  function  in  the  present  case,  as  well  as 
many  other  cases  dealing  with  the  stress  analysis  of  composite 
materials  containing  surfaces,   interfaces,  and  cracks. 

The  solution  of  the  Hilbert  problem  is  obtained  by  using  a 
complex  transform  of  the  type  yL<^  0,5/  where  y  and  q  are  variables 
on  the  real  axis.  It  is  shown  that  this  function  is  an 
eigenfunction  of  the  Hilbert  operator  and  is  orthogonal  for  y, 
varying  in  the  range  between  0  and  ±  »  and  q  between  -oo  to  +  ». 
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PART  I 


Elastic  Green's  Function  for  a  Composite  Solid  Containing 
a  Free  Surface  Normal  to  the  Interface 

by 

V.  K.  Tewary 


1.  Introduction 


In  this  part,  we  describe  the  calculation  of  the  SI  Green's 
function.  We  give  detailed  results  for  the  displacement  and  the 
stress  Green's  function  and  a  discussion  of  their  singularities. 
For  the  purpose  of  illustration,  we  apply  the  method  to  a  simple 
antiplane-strain  problem  for  a  cubic  solid  containing  a  £-5  grain 
boundary  and  also  to  the  corresponding  plane  strain. 

In  calculating  the  SI  Green's  function,  we  have  used  the 
Green's  function  for  an  infinite  bimaterial  composite  as  the 
starting  Green's  function.  A  summary  of  these  functions,  as  taken 
from  two  earlier  papers  [18,19],  has  been  given  in  Appendix  I-A. 
For  calculating  the  SI  Green's  function,  we  need  to  solve  an 
inhomogeneous  generalized  vector  Hilbert  problem.  Our  method  for 
solving  this  Hilbert  problem  has  been  described  in  Appendix  I-B. 
The  evaluation  of  certain  integrals  as  required  for  various 
formulae  in  the  paper  has  been  described  in  Appendix  I-C. 

2.   Calculation  of  the  SI  Green's  Function 

We  consider  a  bimaterial  composite  with  a  planar  interface 
and  a  free  surface  normal  to  the  interface  as  shown  in  figure  1-1. 
This  figure  also  shows  the  coordinate  system  used  in  our 
calculations.  We  take  the  origin  of  our  coordinate  system  at  the 
intersection  of  the  free  surface  and  the  interface.  The  X-axis  is 
taken  along  the  interface  and  the  Y-axis  is  taken  along  the 
surface.  The  Z-axis  is  taken  normal  to  the  plane  of  paper.  We 
assume  that  there  is  no  variation  of  the  field  quantities 
(stresses,  displacement  etc.),  which  enables  us  to  consider  a 
two-dimensional  problem  as  shown  in  figure  1-1.  Although  the  field 
quantities  will  remain  constant  in  the  Z-direction,  they  can  and, 
in  general,  will  have  a  nonvanishing  Z-component. 
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The  notation  in  this  paper  is  generally  the  same  as  in  I  and 
II.  The  solid  in  the  upper  half  plane  (UHP)  and  its  parameters 
will  be  labelled  by  the  superscript  A  and  in  the  lower  half  plane 
(LHP)  by  the  superscript  B.  The  X-,  Y-,  and  Z-  components  of  a 
quantity  will  be  identified  by  subscripts  1,  2,  and  3, 
respectively.  Following  the  usual  notation  in  the  theory  of 
elasticity,  we  shall  denote  the  radius  vector  of  a  point  by  the 
vector  x  with  its  X-,  Y- ,  and  Z-  components  denoted  by  x^  x2 ,  and 
x3 ,  respectively.  The  x2  coordinate  may  also  be  represented  by  y. 

The  Roman  indices  i,j,  k,  etc.  will  take  the  values  1,  2  or  3 
and  will  denote  the  Cartesian  components.  Unless  otherwise  stated, 
we  shall  assume  the  summation  convention  over  repeated  Roman 
indices  but  not  over  repeated  Greek  indices  a,  /3,  A  etc.,  which 
will  also  take  the  values  1,  2  or  3.  We  shall  denote  by 
Greek  i  and  not  by  the  more  usual  symbol  i  which,  as  mentioned 
earlier,  has  been  used  for  labelling  Cartesian  components. 

We  consider  a  case  when  a  line  force  0  is  applied  to  a  solid 
in  the  XY-  plane  at  the  position  vector  x' •  Then  the  displacement 
field  u(x)  at  point  x  is  given  by  the  following  equation  for 
elastic  equilibrium: 

a2  u. 

cikjl   3_    =      0    5(x    _  x<)  6(x    -  x2)    ,  (1) 

Sx^  Sx^ 

where  c  denotes  the  fourth  rank  elastic  constant  tensor  and  5  (x) 
is  Dirac's  delta  function.  We  assume  that  0  is  applied  in  the  UHP 
(x2'  >  0) .  Equation  (1)  is  written  separately  for  UHP  and  LHP  with 
c  and  u  labelled  appropriately  by  A  or  B  for  UHP  and  LHP, 
respectively.   In  the  equation  for  LHP,  the  RHS  will  be  0. 

We  prescribe  the  following  boundary  conditions   (for  i=l-3) : 
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uAi[x1/0]   =  uB.[x1/0]        (0  ^  xx  ^  co),  (2) 


xAi2[Xl/0]  =  TBi2[x1,0]  (0  s  xx  s  »),  (3) 


TAi:L[0,x2]  =0  (0  *  x2  *  co)  ,  (4) 


and 


TBil[0/x2]  =0  (0  *  x2  *  -oo)  (5) 

where  t  is  the  stress  tensor.  The  stress  components  z . .  and  z.0 
can  be  obtained  from  the  displacement  field  by  using 

and 

auk  auk 

Ti2  =  ci21k    ~      +    ci22k    ~       '  (7) 

axl  dx2 

Equations  (2)  and  (3)  are  the  continuity  conditions  which 
imply  that  the  two  solids  are  perfectly  welded  at  the  interface. 
Equations  (4)   and  (5)   are  the  free  surface  conditions  at  x1  =  0. 

Our  object  is  to  solve  eq  (1)  subject  to  the  boundary 
conditions  given  by  eqs  (2)  -  (5)  .  The  Green's  function  will  be 
obtained  from  this  solution  by  taking  the  force  £  equal  to  unity. 
The  corresponding  Green's  function  for  the  case  when  the  unit  line 
force  is  in  LHP  can  be  obtained  by  replacing  x2'  by  and 
appropriate  relabelling  of  the  solids  in  the  UHP  and  the  LHP. 

The  Green's  function,  which  is  a  solution  of  eq  (1)  and  which 
satisfies  the  boundary  conditions  given  by  eqs    (2)    and    (3)  ,  has 
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been  obtained  in  I  (see  also  Appendix  A  of  II)  and  quoted  in 
Appendix  I -A  of  this  paper.  We  shall  use  this  result  to  obtain  the 
solution  of  eq  (1) ,  which  would  also  satisfy  the  boundary 
conditions  given  by  eqs   (4)   and  (5) . 

For  this  purpose  we  shall  use  the  standard  technique  for 
Green's  function  calculations  (see,  for  example,  I  and  II)  and 
apply  a  hypothetical  distribution  of  line  forces  just  outside  the 
region  in  which  the  equation  has  to  be  solved.  We  then  determine 
the  force  distribution  by  requiring  that  the  solution  satisfies 
the  prescribed  boundary  conditions.  In  addition,  we  shall  impose 
the  condition  that  the  total  space  integral  of  the  forces  is  0. 

A 

We   denote   the    force   distributions    by    F  (t)     (0  *  t  *  »)  and 

F  (t)     (-co  <  t  ^  0)    in   the   UHP   and   LHP.    As    shown    in   figure  1-1, 

these    forces    are    applied    just    outside    the    free    surface    at  a 

A 

continuous  set  of  points  in  the  UHP  and  LHP  viz.  F  (t)  at  x=-e 
and  x2  =t  (0  <  t  *  oo )  and  F  (t)  at  x1=-e  and  x2  =  t  (-00  <  t  *  0)  , 
where  e  is  a  small  positive  number  which  is  taken  to  be  zero  in 
the  limit. 

From  the  definition  of  the  Green's  function,  we  can  write  the 
displacement  field  in  UHP  and  LHP;  i.e.,  a  solution  of  eq  (3),  in 
the  following  form: 


o 


+ 


—00 


(8) 


and 


(x;t) FA(t)dt 


0 
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(9) 


—  00 


where  GQ  denotes  the  Green's  function  for  an  infinite  bimaterial 
composite  which  was  derived  in  I.  In  eqs  (8)  and  (9),  we  have  not 
explicitly  written  -e ,  the  x.^  component  of  the  position  vector  of 
the  second  variable  in  the  arguments  of  the  Green's  functions 
under  the  integral  signs,   for  reasons  of  notational  brevity. 

A 

It  may  be  emphasized  that  the  so-called  force  functions  F  (t) 
and  F  (t)  are  arbitrary  functions  which  have  to  be  determined  so 
that  the  boundary  conditions  given  by  eqs  (2) — (5)  are  satisfied. 
These  functions  have  the  dimensions  of  force  but  are  not  physical 
forces  applied  to  the  solid. 

Using  the  result  given  in  eqs  (A.1)-(A.4)  of  Appendix  I -A,  we 
obtain  the  following  expressions  for  the  displacement  field  from 
eqs   (8)   and   (9) : 


a 


a/3 


0 


00 


JT 


1 


0 


Re 


zA+  c-p^  t)     FA(t)  dt 


7 


„0 


+  i 

7T 


Re 


—oo 


ZA,    .    ^III        ,   A  B    .  , 

1  («)  ln(za     e~p/3  t} 

■a/3 


FB(t)  dt 


and 


(for  0  s  x2  s  »)  (10) 


u  (x)  = 


1  n     V"     -B.    .    _II   ,    ,-B     -A  , , 
-  ^  Re  2^    Z  (a)  2/3     ln(za~  p0  ^ 
a/3 


00 


1 

TT 


Re 


Y  IB(a)  Q 
-a/3 


II  i  ,-B,  -A  . . 
_  ln(z  +  e-p0  t) 
/3         v  a        ^/3  ' 


F   (t)  dt 


.0 


1 

TT 


Re 


^  JB(a)   ln(zB+  e-pB  t) 


—00 


a 


,B 


F   (t)  dt 


.0 


+  i 

TT 


Re 


V"  -B,    .    .IV  ,    ,-B^  B 
>     1    (a)  ln(z~+  c_Po  t) 


/3 


—00 


La/3 


F  (t)  dt 


(for  -oo  <  x2  *  0)  ,  (11) 


where  Re  denotes  the  real  part,  overhead  bar  denotes  complex 
conjugate  and  other  symbols  have  been  defined  in  Appendix  I-A. 


Equations  (10)  and  (11)  will  satisfy  eq  (1)  everywhere  in  the 
region  of  solution,  since  the  line  of  application  of  the  forces  is 
just  outside  this  region,  which  is  ensured  by  the  presence  of  c. 
In  the  end  we  shall  take  the  limit  c  =  0. 
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Further,   since  the  Green's  function  used  in  writing  eqs  (10) 

and   (11)    satisfies  the  continuity  conditions  at  the  interface  as 

given  by  eqs    (2)    and    (3)  ,   the  solution  as  given  above  will  also 

A  B 

satisfy  these  conditions  for  all  values  of  F  (t)  and  F  (t)  .  Our 
main  task,  therefore,  is  to  determine  these  two  force  functions  so 
that  the  free-surface  boundary  conditions  as  given  by  eqs  (4)  and 
(5)  are  also  satisfied.  The  resulting  solution  will  then  satisfy 
all  the  prescribed  boundary  conditions. 

In  addition,  we  impose  the  following  condition  on  the  force 
functions  to  ensure  that  the  displacements  are  single  valued  [23] 
everywhere  except  possibly  at  the  cut: 


FA(t)   dt  +  FB(t)   dt  =  0  .  (12) 

0  -co 

As  described   in   I   and   II,    a   simple  method   for  writing  the 
stress  component  x.-   from  eqs   (8)   and   (9)    is  to  differentiate  the 

log  functions  in  these  equations  with  respect  to  their  arguments 

A  B  A  B 

and  replace  the  matrix  y  '     by  the  matrix  Q  '    ,  where  Q  is  defined 

below  (for  superscripts  A  or  B) 


fi.^a)   =  L1^  (a)   rkj(a)  (13) 


and 


L1..  (a)  =  c.1iv  +  p„  c.10V  .  (14) 


The  stress  components  t^2  can  similarly  be  obtained  by 
replacing  y  by  a  where  a  has  been  defined  in  Appendix  I-A.  In  what 
follows  we  shall  denote  the  three  components  t ^  by  the  vector  T^. 
Similarly,  T2  will  denote  the  components  t^2- 
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Thus,    we   obtain  the   following  expressions  for 


(8)   and  (9)  . 


from  eqs 


a 
a/3 


+  L 


I  ,V,  h-(!^L]+5V)„M^L] 
«e  L  pp  1  p£  J 


and 


-  L 


I 

a/3  l 


P/3 


-A 


if 


(15) 


3 


Tj(x)  =  - 


a/3 


+  L 


r  -B  B 

a0   L  p~ 


+  L 


a  Pa 


-B  , 

z  +  c 
a 


10 


-  L 


I 

a/3 


KBB(a,/3)  H 


-B 

Bf  za  ) 


v  P 


B 

/3 


+  KBB(a,/3)  H 


B 


B 


P/3 


(16) 


where 


3A<Z>  =  WZ 


"M  FA(t) 


t  -  z 


dt  , 


0J 


(17) 


5B(Z)  =  2WT 


-°  FB(t) 


t  -  z 


dt  , 


—  00 


(18) 


K^ia.P)  =  QA(a)        /  pA  , 


(19) 


KAB(a/^)  =  QA(oO  Qp11/  p*  # 


(20) 


KBA(a^)   =  fiB(a)   Q*1/  pA  , 


(21) 


KBB(a^)  =  nB(a)  Q*V  /  pB  , 


/3 


(22) 


,A,B,  .  oA,B,  .  /  A, B 
J  '    (a)   =  Q  («)/Pa 


(23) 


and  z,  in  eqs  (17)  and  (18)  ,  is  any  complex  variable.  We  have 
neglected  e  in  the  arguments  of  the  H-f unctions  in  eqs  (15)  and 
(16)  except  in  the  third  term  on  the  RHS,  because  these  are  the 
only  terms  which  will  be  sensitive  to  e  in  the  limit  c=0. 


A  B 

We  now  equate  u.^x)   and  z.  ,  (x)   to  zero  at  x    =  0.  At  x..  =  0, 
A  B 

the  variables     z      and      z      assume  the  following  values: 
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A.B  A.B 
z         =  p        y  , 


(24) 


where  we  have  denoted  the  x^-  component  of  x  by  the  variable  y 
Using  eq  (24)   in  eqs   (15)   and  (16) ,  we  obtain  the  equations 

(i)   For    0  ±  y  2=  w 

TA  HA(y-ie)   +  FA  HA(y+ie) 


a(3  >- 


K^Ca,^)  H^a^y)  +  KAa(a,/3)  HA(5/vfty) 


a/3-1 '  -AA 


a/3- 


"I 


af3 


KAB(a,P)  HB(bogy)  +  KAB(a,p)  HB(b^y) 


=  -^TReZ    ^<a>    ^  Pa/Pa]"'  i 


a 


a0 


(25) 


and 


(ii)   For  -oo  <  y  <  0 


I 

a£  l 


KBA(«,0>   HA(d^y)   +KRa(a//3)  HA(d^y) 


a/3-7 '  -BA 


a/3- 


+  TB  HB(y-Le)   +  TB  HB(y+Le) 
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-z 

af3 


KBB(«,0)  HB(eaey)  +  KBB(a,(3)  HB(eagy) 


'a/3- 


■BB 


-  1 


a0 


(26) 


where 


rA'B  =  £  JA'B(a)  , 
a 


A,  B 


(27) 


P(a,/3)  = 


JA(«)  Qj£  » 


(28) 


R(a,/3)   =  JB(a)  Q*1  , 


(29) 


A    ,  -A 
aa0  =  P<x  /  P/3  ' 


(30) 


b«0  =  Pa  /  P^ 


(31) 


B 


=  Pa  /  P 
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(32) 


-B 


B 


e  _  =  p    /  p_ 
a/3        a  '  ^/3 


(33) 


In  eqs   (25)   and   (26) ,  we  have  used  the  fact  that,    in  view  of 

eq  (A. 20),  the  sign  of  the  imaginary  part  of  the  terms  of  the  type 

e/p^   (for  superscripts  A  or  B)   will  be  negative.   In  the  limit  e  = 

0,  the  functions  H(y  ±  e/Pa)  will  be  sensitive  only  to  the  sign  of 

the  imaginary  part  and  not  its  magnitude,   which  approaches  0  any 

way.    Accordingly,    in   the    argument    of    the    H-f  unctions,    we  have 

written  e/p    as  -lc  and  e/  p    as  +lc. 
' *a  '  *a 


It  may  be  remarked  that  the  quantities  J,  K,  T,  P,   etc.  (with 
various   subscripts   and   superscripts)    are   3x3    square  matrices, 
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whereas    H,     F,     and    <£    are    3x1    column    matrices    or  vectors. 

Equations  (25)   and  (26)   together  represent  a  6  x  6  matrix  equation 

or  a  set  of  6  simultaneous  equations   in  6  unknowns.   The  unknowns 

A  B 

are  the  three  components  each  of  F  (t)  and  F  (t)  .  These  equations 
constitute  a  generalized  inhomogeneous  vector  Hilbert  problem.  Its 
solution  has  been  described  in  Appendix  I-B. 

Following  the  method  described  in  Appendix  I-B,  we  write  the 
two  force  functions  in  terms  of  their  complex  transforms  as  given 
by  eqs  (B.25)  and  (B.26).  In  Appendix  I-B,  the  functions  F(t)  and 
H(z)  are  assumed  to  be  vectors  having  two  components.  In  the 
present  case,   these  are  vectors  having  6  components.   However,  the 

formulae  given  in  Appendix  I-B  are  formally  valid  in  the  present 

ABA  B 
case  as  well,    provided  we  regard  F  (t)  ,    F  (t)  ,    H  (z)  ,    and  H  (z) 

themselves  as  3  component  vectors. 

Thus,  as  in  eqs  (B.27)  and  (B.28),  we  write  the  following 
representation  for  the  H-f unctions: 


TA 


r  (z) 


V^qjz^"0,5  dq  (34) 


—00 


and 


00 


HB(z)  = 


TTB  /  x  tq-0.5  ,  ,.c. 
V  (q)z  M  dq  ,  (35) 


where  (q)  and  V  (q)  are  the  unknown  functions  of  the  real 
variable  variable  q,  which  have  to  be  determined.  For  each  q, 
these  are  3x1  column  matrices  like  the  H-f unctions. 


Once  the  H  functions  are  known,  the  force  functions  can  be 
obtained  in  a  simple  manner  in  terms  of  the  discontinuity  in  the  H 
function  over  the  real  axis  by  using  the  Plemelj  relation  (see, 
for  example,    [20])   as  given  by  eq  (B.42). 
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We  now  substitute  for  the  H  functions  from  eqs  (34)  and  (35) 
in  eqs  (25)  and  (26)  and  use  the  method  described  in  Appendix  I-B 
for  deriving  eqs  (B.31)  and  (B.32).  This  would  transform  each 
element  of  the  matrices  in  eqs  (25)  and  (26)  into  a  function  of  q 
alone.  As  remarked  earlier,  we  can  represent  eqs  (25)  and  (26)  by 
a  single  6x6  matrix  equation.  Accordingly,  we  regard  V(q)  as  a 
single  6x1  column  matrix  or  a  six-dimensional  vector.  Thus,  eqs 
(25)  and  (26)  ,  after  the  transformation  can  be  written  in  the 
following  matrix  representation: 


M(q)   Y(q)   =  N(q)/E(q) , 


(36) 


where  M(q)  and  N(q)  are  respectively  the  6x6  and  6x1  matrices 
formed  by  the  transformed  elements  of  the  LHS  and  RHS  of  eqs  (25) 
and  (26) . 


In  order  to  represent  the  matrices  M(q) ,  N(q)  and  V(q) ,  it 
would  be  convenient  to  introduce  a  block  matrix  notation  in  terms 
of  A  and  B  which  correspond  to  the  solids  in  UHP  and  LHP 
respectively.  According  to  the  block  scheme,  we  first  write  the 
6x1  matrix  V(q)  in  terms  of  two  3x1  blocks  V^q)  and  VB(q)  and 
then  write  M(q)  and  N(q)  in  terms  of  the  same  block 
representation.  Their  matrix  elements  will  thus  be  denoted  in 
terms  of  the  block  matrices  MAB   ,  ,   N^,      etc.   For  notational 

brevity,  we  shall  not  show  their  functional  dependence  on  q  except 
where  explicitly  needed. 


The  block  structure  of  eq  (36)   can  thus  be  written  as 


M 

-AA 


M 

^  -BA 


M 

-AB 


M 

-BB 


1 

>  YB(q), 

E(q) 

-A 


-B 


(37) 


where 
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=  -rA  exp(-27rq)   +  fA 


+ 


Z 

ccfi 


(38) 


a/3>- 


a/3 


(39) 


M 

-BA 


=  Z 


a/3  L 


(40) 


MBB  =  LrBexp(Trq)-LfBexp(-7rq) 


-z 

a/3 


(41) 


N 


A        exp  (-2Trq) 
2n 


a 


JA(a) 


'    Ay  h" 

Pa/Pa 


u 


-A,    %  f -A  /-Al 


+  J"(a) 


/ 


exp(-2Trq)  V" 
271 

a/3 


P  (a,/3) 


r-A/  a)  u 

VP«J 


+  P(a,/3) 


r  A/-ai^ 


<t   ,  (42) 


N 


B    _  l exp  (-Tig) 
27T  E(q) 


z 

a/3 


R  (a,/3) 


(-A  ,-B) 

V     1  J 


+  R(a,/3) 


i   /  (43) 


and 


E(q)   =  1  +  exp(-2Trq), 


l±  =  tq  -  0.5  . 


(44) 
(45) 
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Using  the  same  representation  as  given  in  eq  (37)  ,  the  block 
structure  of  the  matrix  H(z)  is 


H(z)  = 


f  A 
HA(Z) 

I  HB(z)J 


(46) 


Equation  (37)   can  be  solved  for  V(q)   by  matrix  inversion: 


-1 


V(q)   =  M       N  /  E(q)  . 


(47) 


For  later  use,  we  shall  also  write  eq  (47)  in  block  representation 
corresponding  to  eq  (36)  as 


I  YB(q)J 


E(q)D(q) 


c 

1  -BA 


Sab 


Sbb 


> 

'  1 

Sa 

J 

(48) 


where  C  are  the  block  elements  of  the  matrix  of  cofactors  of  M  in 
the  same  block  representation  and  D(q)   is  the  determinant  of  M. 

The  matrix  of  cofactors  is  proportional  to  the  inverse  of  M 

which,  in  the  present  block  representation,  can  be  formally 
expressed  as 

1  StofcrtW*)  =  1  ?Wi>saB<q>  "  D(q)  5ab  -  (49> 

a  a 

where  6,_,   is  Kronecker's  delta  and  A,    a.    and  B  are  block  indices 
AB  '  ' 

which  label  the  blocks  in  the  corresponding  matrix  representation. 
As  in  the  case  of  M  and  N  matrices,  we  shall  show  the  functional 
dependence  of  the  elements  of  C  on  q  only  where  explicitly  needed. 
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In  general  the  inverse  of  M  in  eq  (48)  has  to  be  obtained 
numerically,  and  there  is  no  need  to  calculate  the  cofactor 
matrix.  The  two  H  functions  can  then  be  calculated  from  eqs  (34) 
and  (35)  by  integrating  over  q.  The  integrals  over  q  in  these 
equations  can  be  calculated  numerically  in  a  straightforward 
manner.  It  is  also  possible  to  evaluate  the  q  integral 
analytically  by  using  contour  integration.  For  this  purpose  eq 
(49)  is  more  convenient.  This  would  give  a  series  representation 
of  the  H-f unctions,  which  will  be  derived  in  the  next  section. 

Equation  (47)  gives  the  particular  solution  of  Hilbert's 
equation.  To  obtain  the  general  solution  of  Hilbert's  equation,  as 
described  in  Appendix  I-B,  we  need  to  add  the  following 
homogeneous  solution 

_  • 

HQ(Z)   =  Y  h(Qr)    zLQr"°*5,  (50) 
r 

where  Q    is  a  root  of  the  determinantal  equation 


D(Qr)  =  0,  (51) 

and  h(Qr)  is  an  eigenvector  of  M(q)  corresponding  to  zero 
eigenvalue  for  q  =  Qr;  it  is  a  solution  of  the  matrix  equation: 

M(Qr)h(Qr)  =  0  .  (52) 

In  analogy  with  vectors  H(z)  and  V(q)  the  two  block 
components  of  h(Qr)  will  be  labelled  by  the  superscripts  A  and  B. 
According  to  the  theory  of  homogeneous  equations,  one  of  the 
components  of  h  will  be  arbitrary.  This  provides  a  set  of 
arbitrary  constants  (one  for  each  value  of  Qr)  ,  which  have  to  be 
chosen  so  that  eq  (12)  is  satisfied.  These  constants  can  also  be 
chosen  to  satisfy  any  other  boundary  condition  on  the  stress  field 
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such  as  those  which  may  be  required  at  the  outer  surfaces  in  a 
finite  solid. 


Thus,    finally,   we  can  write  the  H-  functions  which  give  the 
complete  solution  of  eqs   (25)   and  (26)   as  follows: 


H(z) 


HA(z) 
HB(z) J 


(53) 


where  the  subscript  p  denotes  the  particular  solution  which  is 
given  by  eqs  (34)  ,  (35)  and  (47)  .  The  homogeneous  solution, 
labelled  with  the  subscript  0,   is  given  by  eq  (50) . 


Henceforth,  for  notational  brevity,  we  shall  omit  the 
subscript  p  from  the  particular  solution  except  where  we  need  to 
distinguish  it  from  the  homogeneous  solution.  Accordingly,  H(z) 
will  refer  to  the  particular  solution  for  the  H-f unction. 


The  expression  for  the  H-functions  as  given  above  is  the  main 
result  of  this  paper.  The  stress  distribution  in  the  solid  is 
given  by  eqs  (15)  and  (16)  in  terms  of  the  two  H-functions.  The 
displacement  field  is  given  by  eqs  (10)  and  (11)  in  terms  of  the 
logarithmic  integrals  which  can  be  obtained  from  the  H-  functions 
by  a  simple  integration  over  z.  Equations  (10)  and  (11)  give  the 
required  Green's  function  when  £  is  taken  to  be  a  unit  force. 

In  the  next  section,  we  shall  derive  a  series  representation 
for  the  H-functions  and  discuss  their  singularities. 

3.   Series  Representation  of  the  H-Function 

A 

In  the  previous   section  we  obtained  the  two  H   functions  H 
and   H  ,    as    integrals    over   q.    In   this    section   we    shall  obtain 


19 


A  B 

series  representations   for  H     and  H     in  the  Cartesian  space  and 

discuss  their  singularities.   The  singularities  in  the  H-f unctions 

give   the   singularities   in  the   stress.    We   shall   also   obtain  the 

A  B 

space  integral  of  the  force  functions  F  and  F  in  order  to  ensure 
that  eq  (12)   is  satisfied. 

A  B 

First  we  evaluate  H  (z)   and  H  (z)   from  eqs   (34)   and   (35)  for 

A 

a  complex  variable  z.  Using  eq  (48),  we  write  H  (z)   in  the  form 


HA(z)  = 


-co' 


"m  T-.A,   .      iq-0.5  , 

£  <q>  2       dq  (54) 

E(q)  D(q) 


where 


?A(q)   =   [^AA(q)   ^A(q)   +  ^AB(q)   ^B(q)]  (55) 

g 

with  a  similar  expression  for  H  (z) . 

The    q-integral    in    eq    (54)    can    be    evaluated    by   using  the 

method  of  contour  integration  as  described  in  Appendix  I-C  for  eq 

(C.16).  From  the  result  given  by  eq  (C.24)  we  obtain  the  following 

A  B 

expressions  for  H  (z)  and  H  (z)    (in  the  limit  e  =  0) 


HA(z)  =Y^Y}*A{n,e)   2   (n+1)"C  ^(«,n) 


oc  n 


(r,e)  z 


lQ  -0.5+e 
r 
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a£   L  n 


,-(11+1) -e  £  (af,#n) 


r—^  iQ  -0.5+e 

+  )    ?AA(r'c)    Z     r  £2(«/3/LQr) 


a/3   L  n 


,-(n+l) -e 


£3(a/3,n) 


^  lQ  -0.5+e 

+  2^  ^AB(r'C)    Z     r  g3(aP/t-Qr) 
r 


sz  +  Vz> 


(56) 


and 


HB(z)  =  ]T  ^BA(n'e)   z"(n+1)   C  ^(^n) 


L  n 


- ,  lQ  -0.5+e 


a/3   L  n 


,-(n+l)-e 


£2(a/3,n) 
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^— i  lQ  -0.5+c 

+  2^  ^BA(r'C)    2     r  g2(a/3,(.Qr) 


a/3  L 


-(n+l)-c 


£,(*P»n) 


LQ  -0.5+C 

2^  *BB(r,e)   z    r  g3(a/3,LQr) 


sz  +  RH(z)f 


(57) 


where  Re  denotes  the  real  part  and  other  symbols  are  as  defined 
below: 


-AA(n'C)   =  -LexP("2TrLC)   9aa(qn+LC) /D(qTi+2tc)  , 


4n 


n 


-AB(n/C)  =  L(_1)  exP("TrLe)2aR(^+l-e) /D(q„+2Le)  / 


ABv^n 


ln 


$BA(n,e)  =  -Lexp(-27TLC)   CRa  (q^+Lc) /D(qri+2(,c)  , 


BA  ^n 


ln 


n 


$BB(n,c)  =  l   (-1)   exp(-TTLc)   gBB  (qn+Lc )  /D  (qn+2  lc  )  , 


(58) 
(59) 
(60) 
(61) 


2TTLexp[-2Tr(Qr-LC)  ]  9AA(Qr-LC) 
E(Qr-2LC)   D'  (Qr) 


(62) 


*AB(r,c)   =  - 


27T    eXp[-7T(Qr-LC)  ]  CAB(Qr-LC) 

E(Qr-  2  L  c )    D'  (Qr) 


(63) 
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27TL    exp[-27T(Qr-LC)  ]  9BA(Qr-Le) 


E(Qr-2LC)   D' (Qr) 


(64) 


*BB(r,e)  = 


2n  exp[-Tr(Qr-Le)  ]   2gg (Qr~LC) 


E(Qr-2te)   D' (Qr) 


^(ocfn)  =  2  Re  JA(a)    [pA  /  pA] 


n+e 


I^OMQ,.)  =  J*(«)    [p*  /  p*] 


-LQr~0 . 5-e 


-lQ  -0.5-e 
r 


+  SN«>  [p*  /  i£] 

£2(<x/9,n)   =  2  Re  P(a,/3)    [pA  /  pA  J 


n+e 


(65) 


(66) 


(67) 


(68) 


£2(af3,LQr)  =  P(a,/3)    [pA  /pA  ] 


-lQ  -0.5-e 
r 


+!(«,(» [P*  /  p*] 


-lQ  -0.5-e 
r 


(69) 


£3(«0,n)  =  2  Re  R(a,|3)    [pA  /  p*  ] 


n+e 


(70) 


Z3(otfl,LQr)  =  R(a,/3)  [pA  /p*] 


-lQ  -0.5-e 
r 
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(ot,0)  [pA  /  p*] 


-lQ  -0.5-c 
r 


+  R 


(71) 


The   vectors  £2 ,    £3,    £lf    and    £2    are    in   units    of    £/27T.  The 


factor  sz  is  +1  if  the  integration  contour  was  chosen  to  be  in  the 
UHP  (that  is, if  the  expansion  parameter  is  greater  than  1)  and  -1 
otherwise.  In  deriving  the  above  equations,  we  have  taken  the 
following  values  of  the  branches  of  the  complex  exponential 


As  described  in  Appendix  I-C  [see  eqs   (C.22)   and  (C.23)],  the 

allowed  values  of  n  and  Qr  in  each  sum  depend  upon  the  magnitude 

of  the  expansion  variable.    For  example   in  the   first  term  on  the 

RHS  of  eq    (61)  ,    n  will  take  0  or  any  positive   integral  value  if 
A  A 

mod   (z  p  /  pa)    is  larger  than  unity  and  negative  integral  values 

otherwise.  Thus  n  and  also  of  course  Qr  depend  upon  the  indices  of 

the  outer  sum.   The  factor  s     in  eqs    (56)    and    (57)  ,   which  defines 

z 

the  sign  of  each  expansion  term,  also  depends  upon  the  indices  of 
the  outer  sum.  Because  of  this  dependence  we  cannot  write  the 
above  equations  in  a  compact  form  for  a  general  value  of  z. 

However,  in  the  two  limiting  cases  when  mod  (z)  ->  0  or  mod 
(z)  ->  oo,  n  as  well  as  Qr  become  independent  of  a  or  jS.  In  such 
cases  the  above  equations  are  considerably  simplified.  Since  these 
cases  are  of  substantial  interest,  we  give  below  the  matrix 
representations  for  the  two  H  functions  in  the  AB  block 
representation  in  the  two  limiting  cases: 


exp(TTL/2)   =  i  and  exp(-7ri,/2)  =  -l  . 


(72) 


H(z)   =±£  $(n,e)  £(n,e)z 


-(n+1) -e 


n 


lQ  -0.5+e 


±  BH(Z) » 


(73) 


r 
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where  $  and  *  are  square  matrices,  and  H,  R^,  £  and  £  are  column 
vectors  with  elements  labelled  by  A  and  B  in  the  same  block 
representation  as  defined  by  eq  (37).  As  in  eq  (C.24),  the  plus 
sign  is  to  be  taken  if  mod(z)  >  1  and  negative  if  mod(z)  <  1.  The 
matrix  elements  of  $  and  ¥  have  already  been  defined  in  eqs  (58) 
— (65).  The  elements  of  £  and  £  are  defined  as 

^A(n)   =  I  Ma'n)   +  I  £2(«0'n)L  (74) 

£B(n)   =  £    £3(«/3,n)  ,  (75) 
a/3 

^A(LQr}   =  1  £1(«'LQr>  +  1  22(a,3'LQr)'  (76) 
a  a/3 

^B(LQr)   =  I    M0^'^^  '  (77) 
aft 

The  dependence  of  the  various  parameters  in  eqs  (74) -(77)  on 
e  is  not  shown  for  reasons  of  notational  brevity.  As  described  in 
Appendix  I-C,  e  can  be  set  to  zero  immediately  for  those  values  of 
n  and  r  for  which  there  are  no  higher  order  poles.  When  the  higher 
order  poles  are  present,  their  contribution  has  to  be  calculated 
separately  by  evaluating  the  limit  at  e  =  0.  In  eq  (73),  this 
contribution  has  been  denoted  by  R„(z) .  Equations  (73) ,  (56)  and 
( 57 )  ,  therefore,  provide  a  convenient  representation  for  the  H 
functions  which  is  formally  valid  whether  higher  order  poles  are 
present  or  not.  In  section  5  we  shall  consider  the 
antiplane-strain  case  in  which  all  qn's  are  also  roots  of  eq  (51) 
and  are  therefore  second-order  poles. 

As  explained  in  Appendix  I-C  and  given  by  eq  (C.22)  and 
(C.23),  the  position  and  the  contribution  of  poles  depend  upon 
whether  the  contour  of  integration  has  been  chosen  in  the  UHP  or 
LHP.    This,     in    turn,     depends    upon    the    magnitude    of    z    or  the 
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expansion  variable.  In  sections  3B  and  3C,  respectively,  we  have 
calculated  the  contribution  of  the  poles  in  the  two  limiting  cases 
when  mod(z)  ->  »  and  mod(z)  -»  0. 

We  see  from  eq  (73)  ,  and  also  from  eqs  (56)  and  (57)  ,  that 
the  H-f unctions  contain  basically  two  series  representations,  one 
in  qn  and  the  other  in  Q  ,  and  a  remainder  term  which  is 
nonvanishing  when  there  are  higher  order  poles  in  the  integrand. 
The  terms  containing  qn  are  actually  those  which  give  the 
particular  solution  of  eqs  (25)  and  (26) .  The  terms  containing  Qr 
are  the  solution  of  the  homogeneous  equation  because,  when 
substituted  in  eqs  (25)  and  (26)  ,  they  lead  to  the  homogeneous 
matrix  equation  (36)    (unless  qn  =  Qr  for  some  n  and  r) . 

Thus,  although  eqs  (56)  and  (57)  are  particular  solutions  of 
the  Hilbert's  equation  in  the  q  representation,  they  may  still 
contain  a  homogeneous  part  which  separates  out  when  integrated 
over  q.  We  shall,  however,  refer  to  these  terms  also  as  part  of 
the  particular  solution  to  distinguish  it  from  the  homogeneous 
part  which  has  to  be  separately  introduced  in  order  to  satisfy  eq 
(12) .  This  part  will  be  discussed  a  little  later  in  section  3D. 

3A.  Behavior  of  H^(z)   for  large  values  of  mod(z) 

We  shall  now  discuss  the  behavior  of  H(z)  given  by  eq  (73) 
for  large  mod(z)  .  In  this  limit,  as  given  by  eq  (C.22),  only  the 
poles  in  the  UHP  will  contribute  to  the  integral.  The 
corresponding  allowed  values  of  n  and  Qr  are  then  given  by  eq 
(C.22).  For  a  moment  if  we  ignore  any  higher  order  poles,  we  find 
from  eqs  (73)  and  (C.22),  that  H(z)  will  consist  of  only  negative 
powers  of  z.  This  shows  that  the  stress  will  vanish  as  mod(z)  ->  «. 

The  first  series  in  eq  (73)  contains  negative  integral 
powers,    whereas    the    second    series,     in    general,     will  contain 


26 


negative  fractional  powers  of  z.  Let  us  write  Qr  in  terms  of  real 
and  imaginary  parts  as 


Q  =  g  +  t.k  , 
*r      ^r  r' 


(78) 


with  the  constraint 


kr  =  Im(Qr)   >  0. 


(79) 


The  second  series  in  eq  (74)   can  then  be  written  in  the  form 


The  subscripts  p2  on  the  vectors  V  and  H  identify  them  as 
belonging  to  the  second  series  of  the  particular  solution.  In  eqs 
(80)  and  (81)  ,  we  have  put  c  =  0  consistent  with  our  hypothesis 
that  the  particular  value  of  Q    is  not  a  higher  order  pole. 

We    notice    from    eq     (81)     that    each  term    falls    off  with 

increasing  mod(z) .    The   exponential   factor  for   each  term  will  be 

nonunity  whenever  gr,  the  real  part  of  Q  corresponding  to  that 
term  is  nonzero. 

This  discussion  is  valid  for  those  terms  which  arise  from 
only  simple  poles.  However,  as  shown  below,  at  least  one  pole, 
that  at  n=0,  is  a  second-order  pole.  For  this  pole,  therefore,  we 
have  to  evaluate  the  proper  limit  at  c  =  0. 

It  can  be  verified  that  D(q)    is  zero  at  q  =  q     where  q  has 


(80) 


[for  mod(z)  -»  oo] , 


where 


Y  2(Qr,c)  =    *(Qr)  £(LQr). 


(81) 
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been  defined  in  eq  (C.25).  For  this  purpose,  we  note  from  eqs 
(38)-(41)/  (84)  and  (85),  and  by  using  the  definition  of  various 
parameters  as  given  in  section  2  and  Appendix  I-A,  that  at  q  =  qQ, 
the  matrix  elements  of  M(q)   and  £(q)   reduce  to 


*AA<V   "  2  Re[  ^  +  CA  2s  ]    •  (82) 


Sm<V  =    "2  Re[  ?  25U  ]  •  <83> 


5ba<V  -   "2  Re[  CB  2"  ]  - 


(84) 


?WV   -     "2  Re  "  CB  2^  ]    .  (85) 


£a<°>=  '  (86> 

£B(0)  "    "  ^ba'^q)   •  (8?) 

The  quantities  Q1  -  QIV  in  eqs  (82) -(85)  have  been  defined  in  eqs 
(A. 7) -(A. 10)  in  Appendix  I-A  .  They  should  not  be  confused  with  Qr 
which  denotes  a  particular  value  of  q. 

Using  eq  (A. 25)  of  II,  we  find  that 


-  ?WV  =  ?WV  -  !?bb<V  -  °<  <88> 

which  shows  that 

D(qQ)   =  0.  (89) 

Since  E(q)  also  has  a  zero  at  qQ  [see  eq  (C.17)],  we  see  that  the 
integrand  in  eq  (54)  has  a  second-order  pole  at  q    =  Q.  =  i/2. 
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The  contribution  of  this  second  order  pole  has  been  evaluated 
in  Appendix  I-C.  Since  this  pole  is  in  UHP,  it  does  not  contribute 
to  the  integral  in  the  limit  mod(z)  -»  0.  For  large  values  of 
mod(z) ,  we  calculate  its  contribution  by  using  eq  (C.29)  as 
described  below. 

By  comparing  eqs  (C.16)  and  (54),  we  note  that  in  the 
calculation  of  (z)  ,  the  function  P(q)  as  defined  by  eq  (55)  for 
q  =  qQ  is  given  by 

EAC3o»  -  [£aa<«o>  SA(q0)  +  £AB<V  HB<V]-  (90> 

In  order  to  achieve  some  simplification,  it  would  be  useful 
at  this  point  to  derive  certain  relations  between  the  block 
elements  of  the  matrix  C(q)  at  q  =  qQ,  From  eqs  (49)  and 
(86)-(87),  we  see  that 

Wo'  £a<°>  "  sABC30)  £B(°>  =  D<V  -  0  <91> 

and 

SbA^0»   £a<0'   "     SBB(q0)   CB(0)   =     0.  (92) 

We  can  derive  the  following  relations  by  taking  the  product 
of  M(q)   and  C(q)   in  that  order  and  by  using  eqs   (49)   and  (88) : 

Saa^o'  +  SBA<V  =  0  <93> 

and 

SAB(q0)    +  £BB(q0)   =  0.  (94) 

Now  we  come  back  to  eq  (90) .  From  the  definition  of  N(q)  as 
given  by  eqs  (42)   and  (4  3) ,  we  note  that 

W  =  "  £a(0)  (95) 
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and 

SB(q0)   =       ^B(0)*  (96> 
From  eqs  (91) ,    (95)   and  (96)  we  find  that  eq  (90)  gives 

PA(qQ)   =  0.  (97) 

Similarly,  we  can  show  that 

PB(qQ)   =  0.  (98) 

Hence  we  need  to  keep  only  the  first  term  on  the  RHS  of  eq  (C.29). 
A  similar  argument  applies  to  the  calculation  of  H  (z)  .  Thus,  in 
the  limit  c  =  0,  we  get  the  result  [for  mod  (z)  -»  «]  : 

R^(Z)   =         L   Z  1  PA/ (qQ)  (99) 


and 


D' (qQ) 


l    z"1  ..B' 


R„(z)   =   —        P     (qn)  ,  (100) 

where  the  primes  denote  differentiation  with  respect  to  the 
argument . 


Equations  (99)  and  (100)  give  the  contribution  of  the 
second-order  pole  at  qQ  to  the  H  functions,  if  there  are  no  other 
higher  order  poles  in  the  integrand  of  the  H-f unctions.  In  case 
some  other  values  of  q  and  Q  coincide,  we  need  to  take  the  limit 
c  =  0  of  those  terms  in  a  similar  manner,  and  they  would  be 
included  in  R  (z)  .  In  those  terms  in  eqs  (56)  and  (57)  ,  in  which 
qn  is  not  equal  to  Qr,  c  can  be  put  equal  to  zero  immediately. 

The  pole  at  q^  as  discussed  above  is  effectively  a  first 
order  pole  because  P(qn)    is  also  zero  along  with  D(qn)    at  q  =  qn. 
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This  is  the  reason  why  logarithmic  terms  cancel  out  in  the  final 
result  given  by  eqs  (99)   and  (100) . 

In  deriving  eqs    (99)    and    (100)  ,   we  have  only  considered  the 

degeneracy  of  the  poles  between  set  I  and  set  II  that  is  between 

qn  and  Qr  as  defined  by  eqs   (C.17)   and   (C.18).   It  is  obvious  from 

eq  (C.17)   that  there  cannot  be  a  degeneracy  among  the  poles  of  set 

I  because  no  two  q  's  can  be  equal. 
.  n 

The  poles  of  set  II  as  defined  by  eq  (C.18)  might  be 
degenerate  among  themselves.  Since  the  matrix  M(q)  is  6x6,  it 
would  have  six  eigenvalues.  The  determinantal  equation  (C.18)  can, 
therefore,  have  up  to  6-fold  degenerate  roots.  If  tj  is  the 
degeneracy  of  a  particular  root,  that  is  if  tj  number  of  Q  are 
equal,  then,  as  given  by  the  Cauchy  integral  theorem,  their 
contribution   to   the   integral  *  will   contain   a   factor   of   the  type 

77-I 

[ln(z)]  .  The  function  H(z)  will  not  diverge  at  00  because  the 
logarithmic  term  will  be  multiplied  by  a  negative  power  of  z. 

3B.  Behavior  of  H  (z)   for  low  values  of  mod(z) 

P 

In  this  subsection  we  shall  discuss  the  form  of  the 
particular  solution  for  H(z)  in  the  limit  mod(z)  ->  0.  The  behavior 
of  H(z)  for  low  values  of  mod(z)  can  be  analyzed  following  the 
same  steps  as  in  the  preceding  subsection  for  large  mod(z).  The 
singularities  near  mod(z)  ~  0  in  the  total  H(z) ,  including  the 
homogeneous  solution,  will  be  discussed  separately  in  section  3E. 

The  expression  given  in  eq  (73)  for  H(z)  is  also  valid  in  the 
limit  mod(z)  ->  0,  but  the  allowed  values  of  n  and  Qr  will  be  those 
given  by  eq  (C.23).  The  first  series  in  eq  (73)  will  contain  only 
positive  integral  powers  of  z.  These  terms  will  vanish  in  the 
limit  mod(z)  ->  0.  The  second  series  may  contain  also  fractional 
powers.  Except  for  those  values  of  Q    which  have  an  imaginary  part 
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greater  than  -0.5,  the  powers  of  z  in  the  second  series  will  also 
be  positive.  This  shows  that,  apart  from  the  exception  given 
above,  the  terms  in  both  the  series  in  eq  (73)  will  approach  0  as 
mod(z)  goes  to  0. 

In  the  second  series,  we  again  write  Q     in  terms  of  its  real 
and  imaginary  parts  as  given  below  [c.f.   eq  (78)] 

Q  =  g  -  Lk  (101) 
*r      ^rp         rp  v  ' 

with  the  constraint 


k      =  -Im(Q  )   >  0.  (102) 
rp  v  r'  v  ' 

The  subscript  p  on  g  and  kr  is  introduced  to  identify  this  term 
as  coming  from  the  particular  solution.  The  resulting  expression 
for  H  ~(z) ,  which  is  analogous  to  that  given  by  eq  (80)  is 

P<£ 

(k  -0.5) 

Hp2(z)  =  1  Yp2(Qr)    [z]     ^  exp[Lgrpln(z)]  (103) 


for  mod(z)  ->  0. 


Now  we  consider  the  contribution  of  the  higher  order  poles 
following  the  same  procedure  as  used  in  the  preceding  subsection. 
It  may  be  verified  that  D(q)  is  zero  at  q  =  Q  =  -i/2.  Since  E(q) 
is  also  zero  at  q  =  q_1  =  -i/2,  the  integrand  has  at  least  one 
second-order  pole  at  q  =  Q  =  -l/2.  Unlike  the  case  for  mod(z) 
->  oo  as  discussed  in  the  previous  section,  P(q_1)  is  not  zero  in 
this  case.  The  contribution  of  this  second-order  pole  which  is 
represented  by  the  remainder  term,  has  been  evaluated  in  Appendix 
I-C  [see  eq  (C.30)]  and  is  given  below: 
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D' (q.x)  L 


P'(q_1)   +  LP(q_1)  Ln(z) 


-0.5  P(q-:L)D"  (q_1)/D'  (q_1)    +  7TP(q_1) 


(104) 


Equation  (104)  shows  the  ln(z)  behavior  of  H(z)  at  low  values 
of  mod(z)  .  In  deriving  this  equation,  we  have  assumed  that  the 
root  Q_1  of  the  determinantal  equation  (C.18)  is  not  degenerate. 
As  remarked  earlier,  if  eq  (C.18)  has  a  7}-fold  degenerate  root  at 
Qr7J/  H(z)  will  have  terms  of  the  type, 


RR(z)  *  [z] 


(k  -0.5) 
v  r7]  ' 


T)  ™  1 

[ln(z)]7?~  exp[Lgr7)ln(z)  ]  , 


(105) 


where  gr^  and  kr^  are,  respectively,  the  real  and  imaginary  parts 
of  Q  .  If  k  is  half-integer  and  g  =0,  then  E(q)  will  also  be 
0  at  q  =  Q  .  This  is  similar  to  the  case  previously  considered 
for  which  Qri^  =  -l/2.  This  will  lead  to  an  additional  power  of 
ln(z)   as  in  eq  (104) . 


3C.  Contribution  of  the  homogeneous  solution 


We  now  consider  whether  the  homogeneous  solution  HQ(z)  has  to 
be  included  in  the  total  solution  or  not.  For  this  purpose  we  need 
to  calculate  the  integral  of  the  force  functions  with  the 
objective  of  ensuring  that  eq  (12)  is  satisfied.  The  integral 
required  for  calculating  the  force  function  from  the  H  function 
has  been  evaluated  in  Appendix  I-C.  However,  we  shall  not  need  the 
explicit  expressions  for  the  force  functions  for  calculating  the 
Green's  functions;  a  knowledge  of  the  H-f unctions  will  suffice.  In 
this  section,  therefore,  we  need  only  to  consider  the  space 
integral  of  the  force  function,  which  also  has  been  evaluated  in 
Appendix  I-C  by  using  contour  integration. 
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As  shown  in  Appendix  I-C,  the  only  poles  which  contribute  to 
the  force  integral  are  those  at  Qr,  which  are  the  zeroes  of  D(q) 
and  which  satisfy  the  constraint  given  by  eq  (C.35).  Each  value  of 
Qr  corresponds  to  a  force  term  containing  LQr~0.5  as  the  power  of 
y  and  a  stress  term  containing  the  same  power  of  y.  Each  of  these 
terms  will  be  the  homogeneous  solution  because,  by  definition, 
D(Qr)  will  be  0  for  all  of  them.  We,  therefore,  choose  additional 
force  terms  such  that  their  integral  is  negative  of  that  given  by 
eq  (C.34)  so  that  the  total  integral  of  the  force  function  is  0  as 
required  by  eq  (12) .  The  corresponding  H-f unctions  which  denote 
the  homogeneous  solution  are  written  in  the  form, 

h£(z)  =  -    £        v£  (Qr)   zLQr"  °-5  ,  (106) 
Qr*<-/2 


5o<z>  =  "   I      vo  <V  zLQr 


lQ  -  0.5 


Qr*L/2 


and 


(107) 


exp  (ttQ  )  PB(Qr) 

V°(Q  )   =     27T  —  ,  (109) 

E(Qr)   D' (Qr) 

with  the  constraint 

0  ^  Im(Qr)   *  0.5  .  (110) 
3D.  Singularities  in  H(z)   at  mod(z)  =  0 

Now  we  discuss  the  singularities  in  H(z)   as  mod(z)  approaches 
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0.  As  given  by  eq  (110),  H^z)  contains  two  main  terms,  one  is  the 
particular  solution  which  is  an  integral  over  q,  and  the  other  is 
the  homogeneous  solution  which  is  a  sum  over  Q  ,  where  Qr  is  a 
root  of  eq  (51)  with  the  constraint  given  by  eq  (110) . 

First,  we  consider  the  homogeneous  solution.  As  described  in 
the  preceding  section,  this  term  will  arise  only  if  eq  (51)  has  a 
solution  QrQ  which  satisfies  the  constraint  given  by  eq  (109) .  The 
additional  subscript  0  has  been  introduced  to  identify  this  term 
as  arising  from  the  homogeneous  solution.  This  singularity  can  be 
expressed  as 

-(k  +0.5) 

H(z)  *  (z)       ru  exp[(,gr0ln(z)  ]  ,  (111) 

where  the  weight  (coefficient)  of  each  singular  term  is  given  by 
eqs  (108)  and  (109)  ,  and  grQ  and  krQ  are,  respectively,  the  real 
and  imaginary  parts  of  Q  n, 


Qr0  =  gr0  +  «.kr0  ,  (ii2) 

where 

0  s  krQ  s  0.5  .  (113) 

We  see  from  eq  (111)  that  the  exponent  of  the  singularity 
varies  between  0  and  -1.  In  addition,  if  the  real  part  of  QrQ  (in 
addition  to  the  constraint  on  its  imaginary  part)  is  nonzero,  H(z) 
will  also  have  a  logarithmic — cos  or  sin[g  ln(z)] — type 
oscillatory  component,  which  is  represented  by  the  exponential 
factor  in  eq  (111) .  This  oscillatory  component  is  of  the  same  type 
which  arises  in  the  stress  and  the  displacement  field  near  an 
interfacial  crack  in  bimaterial  composites  (see,  for  example,  II) . 
Depending  upon  the  material  parameters  of  the  two  constituents  of 
the  composite,   this  term  can  be  significant.   The  singular  term  as 


35 


given  by  eq  (111)  will  of  course  be  zero  if  a  Qr  does  not  exist 
which  satisfies  the  constraint  given  by  eq  (113) . 

Next,  let  us  identify  the  singular  terms  in  the  particular 
solution  for  H(z)  ,  which  is  given  by  an  integral  over  q  in  eq 
(54)  .  As  shown  in  section  3A,  for  mod  (z)  «  0,  this  integral  can 
be  represented  in  terms  of  two  series  and  a  remainder  term  as 
given  in  eq  (73)  .  As  explained  in  section  3B,  the  first  series  in 
eq  (73)  has  no  singularity.  The  second  term  in  eq  (73)  is  a  series 
over  Qr,  where  Qr  are  solutions  of  eq  (51)  with  the  constraint 
that  Im  (Qr)  is  negative  since  these  terms  arise  from  the  poles  in 
the  LHP.  This  series  is  given  by  eq  (103)  . 

We  see  from  eq  (103)  that  only  those  terms  will  be  singular 
at  mod(z)  =  0  for  which  k  <  0.5.  If  the  real  part  of  Q  is 
nonzero,  the  logarithmic  oscillatory  factor  will  also  be  present. 
Thus  we  obtain  the  following  condition  for  a  term  in  the  second 
series  in  eq  (73)   to  be  singular: 


0  ^  krp  s  0.5  .  (114) 

The  corresponding  singular  terms  in  H(z)   are  of  the  type, 

(k  -0.5) 

Hp2(z)  *  (z)     rp  exp[Lgrpln(z) ] ,  (115) 

where  k  has  to  satisfy  the  constraint  given  by  eq  (114) ,  and  the 
additional  subscript  2  in  H  identifies  the  second  series  in  eq 
(73)  as  the  origin  of  this  singularity.  The  exponent  of  this 
singularity  varies  between  0  and  -0.5. 


Finally,  we  consider  the  singularities  in  the  remainder  term 
R„(z)   in  eq  (73),  which  arise  if  there  are  higher  order  poles.  The 

hi 

contribution  of  one  second-order  pole  at  q  =  -i/2,  at  which  both 
E(q)     and    D(q)     are    0,     is    given    by    eq     (104)  .     This    gives  a 
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singularity  in  H(z)   which  varies  as  ln(z).   This  singular  term  has 


(51)  has  degenerate  roots,  H(z)  will  have  higher  order  logarithmic 
singularities  as  given  by  eq  (105) .  In  general,  the  logarithmic 
oscillatory  factor  will  also  be  present  as  shown  in  eq  (105) . 

The  singularities  in  the  H-function  arise  from  the 
homogeneous  part  of  eq  (36)  and  are,  in  general,  independent  of 
the  inhomogeneity ,  that  is  the  RHS  in  eq  (36) . 

4.  Expressions  for  the  Green's  Functions 

In  this  section  we  shall  give  the  final  expressions  for  the 
displacement  Green's  function  and  the  stress-2  Green's  function. 
The  former  gives  the  displacement  and  the  latter  the  stress 
components  t^2  for  a  unit  line  force.  The  stress-1  Green's 
function,  which  gives  the  stress  components  t-1  for  a  unit  line 
force  is  given  by  eqs  (15)  and  (16)  .  We  shall  give  the  results 
only  for  the  case  when  the  unit  line  force  is  in  the  UHP  at  the 
point  (x',y')  with  p    defined  as 


The  corresponding  expressions  for  the  case  when  the  unit  line 
force  is  in  the  LHP  can  be  obtained  by  a  simple  change  of 
variables,  y'  and  p. 

4A.  Displacement  Green's  Function 


no  oscillatory  factor  since  the  real  part  of  Q       is   zero.    If  eq 


P 


a 


=  x'  +  p    y'  . 


(116) 


a 


a/3 
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+  L 


a  u 


IA(a) 


A , 


-A 

f     Z     +    C  N 

-«  


a 


a 


+  L 


-A 

z 


-  L 


-A 

z 


a/3   L  p/3  P/3 


(117) 


(for  0  s  x2  s  oo) 


and 


GB(x,x')  =  -  |  Re  ^    iB(a)   Q*1  ln(z*-  pA  )g 

a(3 


+  L 


a/3 


-B 

z 


EBA(a,/3)   l^f  +  ERa(a,/3) 


A 

1  pe 


;BA 


B 

Z  ^ 
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+  L 
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B ,  .  TTB 
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B^ 
f  z  +  c 
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B 


+  rB(«)  yB 
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V  C 
-B 


-  L 


I 

a/3  <- 


3 


-B 

f      Z  X 


ERR(a,/3)   U"  +  ERR(a,0)  U 


B 


B 

Za 


BB 


(118) 


(for  -oo  s  x2  s  0)  , 


where 


(119) 
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EAB(a,/3)   =  yA(a)  Q1 


/3 


(120) 


EBA(<*,/3)   =  iB(a)   Qj1  , 


(121) 


(122) 


^(z)   =  - 


27TL 


00 


F  (t)   ln(z  -  t)dt 


0J 


dz 


(123) 


and 


yB(z)  = 


2TTL 


FB(t)    ln(z  -  t)dt 


-oo  J 


HB(z)dz  . 


(124) 


In  deriving  eqs  (117)  and  (118)  from  eqs  (10)  and  (11)  ,  we 
have  neglected  terms  which  contain  factors  of  the  type  ln(pa)  if 
they  represent  only  rigid  body  displacements.  The  integrals  in  eqs 
(123)  and  (124)  are  indefinite  integrals.  These  integrals  can  be 
evaluated  by  using  either  the  integral  or  the  series 
representation  for  the  H-functions  which  have  been  derived  in  the 
previous  sections.  It  can  be  verified  by  using  eq  (6)  that  eqs 
(117)   and  (118)  yield  eqs  (15)   and  (16)   for  stress  components  t^^- 


4B  Stress-2  Green's  Function 


The  stress-2  Green's  function  which  can  be  derived  from  eqs 
(117)   and  (118)   by  following  the  procedure  given  in  section  2  or  I 
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and  II,  are  given  below: 


T2(x,x<)  =  -  i  Re  £    aA«x)    [zA-  pA  ] 


-1 


a 


1  _  A/N.I      [A  -Al 

-  -  Re  ^    a  (a)  \za-  PfB  j 

a/3 


-1 


+  L 


(X 


,  za+  c  ^ 


a 


-A, 

z  +  e  n 


+ 


f-A,  .  ,  -A\  UA(  *aT  fc  i 
{2  («>/  Pa}  «  I  -=s-j 


a 


+  t 


-A 

z 


£     Baa^.P)  hA(  -£")  +  5AA(«^)  hA( 

a/3    L  P^ 


-  L 


a/3  L 
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DAB(a,^)  H 
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and 
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B 


B 

z 

a 


(126) 


where 


D^fof, 0)   =  aA(a)         /  p 


-A 


/3 


(127) 


DAB(a,0)  =  aA(ot)  Ql11/ 


s0 


<3 


(128) 


DBA(«^) 


-B .    .    _II  /  -A 

e  («)  / 


(129) 


and 


DBB(<x,/3)  = 


-B.  .  ^IV  /  B 
2  (a)  2/3    /  P/3 


(130) 


In  the  next  two  sections,  we  shall  illustrate  the  formalism 
developed  in  the  preceding  sections  by  applying  it  to  two  simple 
cases:  the  antiplane-strain  problem  and  the  plane  strain  problem 
in  a  cubic  solid  containing  a  Z-5  tilt  grain  boundary. 
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5.   Application  to  the  Antiplane-Strain  Problem 


For  illustration,  we  apply  the  formalism  developed  in  this 
paper  to  the  simple  case  of  antiplane-strain  in  a  bimaterial 
composite  containing  a  free  surface.  This  case  is  of  little 
practical  interest  but  has  the  important  advantage  of  simplicity. 
It  can,  therefore,  serve  the  purpose  of  illustrating  the 
mathematical  technique  for  calculating  the  SI  Green's  function 
without  unnecessarily  clouding  the  essential  features  of  the 
technique.  The  experience  gained  in  treating  this  simple  case 
should  be  useful  to  readers  in  solving  a  realistic  problem.  The 
calculation  of  the  SI  Green's  function  for  a  more  realistic  case 
has  been  given  in  the  next  section. 

In  the  antiplane-strain  case,  all  the  3x3  matrices  such  as  y, 
a,  Q  and  of  course  the  Green's  function  itself  reduce  to  lxl 
matrices,  that  is,  pure  numbers.  Further,  in  this  case  the  indices 
a,  |3,  etc.,  will  take  only  one  value,  that  is,  a  =  j3  =1.  We  assume 
the  solids  A  and  B  which  constitute  our  composite  model  to  be 
cubic.  This  is  the  same  model  which  has  been  considered  in  I  and 
II.  Values  of  various  parameters  for  this  model  are  given  below. 
Some  of  these  have  been  calculated  in  I  and  II;  others  can  be 
obtained  quite  simply  by  following  the  same  method. 

P«  =  P/3  =  L    '  (131) 
rA'B(«)   =  *g'B  =  l/2dAB   ,  (132) 

crA'B(a)  =  (AB  =  i/2        ,  (133) 

Q1  =  Q1  =  QIV  =     QIV  =  d     ,  (134) 
a        s        a  s  '  v  ' 

Q11  =    Q11  =  2dD     ,  (135) 
a  s  B  v  ' 

Q111  =     Q111  =  ~2dA    ,  (136) 
a  s  A  v  ' 

nA,B  =  1/2    ;  (137) 
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KAA(tt/<3)    =  Ld/2;     KAB(a,<S/    =  LdA  ' 
KBA(«,^)   =  R(a,/3)   =  Ldg  , 

KBB(a,^)  =  P(a,/3)  =  -0.5  Ld  , 
jA,B  =  r  =  _L/2  ^ 

EAA  =  d/2dA  '   EAB  =  -1  ' 

EBA  =  1       EBB  =  d/2dB  ' 
aa/3  =  -1  =  exp(TTL)    ;  =  -1  =  exp(-Tri)    ,  (138) 

ba/3  =  da(3  =  1  ' 
d  =  dA  -  dB     ,  (139) 

dA  _  =  cA,B  /  [  cA    +  cB    ]    ,  (140) 

A,B  44      '      L        44  44    J     *  v  ' 

dA  +  dB  =  1   .  (141) 


Substituting  the  values  of  the  constants  in  eqs  (131) -(141) 
into  eqs  (117)  and  (118) ,  we  obtain  the  following  expressions  for 
the  displacement  Green's  function  for  a  unit  line  force  applied  at 
x'   in  the  UHP  : 

i.  x  in  UHP 

G(x,x')   =  (-l/2TTdA)   Re[ln(z-p)   +  d  ln(z-p)] 


(142) 
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ii.  x  in  LHP 


GB(x,x')  =  (-1/7T)  Re[ln(z-p)]  +  t-[uA[-  \  )   +  U*(  f  ]] 


-  S^l  !)♦"*(-? )]  • 


(143) 


Following  the  procedure  given  in  section  2,  or  directly  from 
eqs  (15)  and  (16)  ,  we  obtain  the  following  expressions  for  the 
stress  (r   )  Green's  functions: 


TA(x,x')   =  -  Re  — 

x  2  71  L  (z 


+ 


(z  -  p)  (z  -  p  ) 


-|—   [  HA[(z+e)/i]   -  HA[-(z  +  C)/L 


[  HA[-z/L]   -  HA[z/L]]    +       dA   [  HB[z/L]   -  HB[-z/c] 


and 


TB(x,x')  =  -  Re   [      1     _^  ]     -  dB  [  HA[-z/L]   -  HA[z/i] 


|  [  HB[-(z+e)/i]  -  HB[(z+e)/t 


EB[z/l]   -  HB[-z/l]  j    ,  (144) 
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where,   in  view  of  eqs  (131)   and  (A. 5), 


z  =  x    +  l  x     =    x  +  lv. 
12  J 


(145) 


As  in  section  2,  we  obtain  the  following  equations  for  the 
two  H-f unctions  by  using  the  boundary  condition  given  by  eqs  (4) 
and  (5)   or  directly  from  eqs  (25)   and  (26) : 


HA(y  +  lc)   -  HA(y  -  lc)  j  = 


27TL 


k  y  +  Lp 


y  -  Lp 


and 


+  d 


L  y  +   tp  y  -   Lp  JJ 


(146) 


[  HB(y  +  lc)   -  HB(y  -  lc) 


Ly-Lp  y+LpJ 


(147) 


Following  the  method  given  in  Appendix  I-B  and  the  steps 
given  in  section  2,  we  obtain  the  following  solution  of  eqs  (146) 
and  (147)  . 

V^q)   =  s(-Lp)   -  s(Lp)   +  d  [  s(-Lp)   -  S(LP)    ],  (148) 

VB(q)   =  2  dBL  exp   (irq)    [  s(-Lp)    -  s(lP)    ]    ,  (149) 
HA(z)   =     f  Aq)    zLq"0*5  dq   ,  (150) 

—oo 

and 
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H*(Z)    =     f  VB(q)    zLq-°'5  dq   ,  (151) 


where 


s(p)   =  2jf-  [exp(-27rq)  /  E2(q)j   p  Lq  °*5   .  (152) 

Now  we  consider  evaluation  of  the  integral   in  eq    (150)  with 
(q)   given  by  eq   (148) .   The  integrand  has  second-order  poles  at 
all  values  of  q    [see  eq   (C.17)   for  definition  of  q  ].  We  can  use 

n  n 

eq   (C.24)    for  evaluating  this  integral  by  putting  all  Q    =  q  and 

r  n 

accounting  for  the  second-order  poles  as  described  in  Appendix 
I-C.  This  would  give  a  series  in  q    which  can  be  easily  summed.  We 

n 

can  also  evaluate  this  integral  directly  by  making  the 
substitution  t  =  exp(-2irq).  The  limits  of  integration  over  t 
become  0  to  »,  The  integral  can  then  be  carried  out  by  using  the 
same  contour  as  used  in  Appendix  I-C  for  evaluating  I(z)  in  eq 
(C.l).  The  integral  for  H  (z)  in  eq  (151)  can  be  evaluated  in  a 
similar  manner.  Thus  we  obtain  the  result 


Hp(z)   =  ~  2nT  [  S^'"^)   "  S(Z'LP) 


+  d       S(z,-lP)    -  S(z,Lp)  (153) 


and 


HP(Z)    =  "  WZ   [S<Z'LP)    ~  S(z,-Lp)J    ,  (154) 


where 


S(z,p)  =  2ir(z  -  p)ln  (Z/P)    *  (155) 


46 


In  eqs  (153)  and  (154)  ,  the  subscript  p  on  the  H-f unctions 
indicates  that  they  are  particular  solutions  of  the  Hilbert's 
equation.  We  have  to  add  the  homogeneous  solution  to  the 
particular  solution  unless  the  corresponding  force  functions 
satisfy  the  condition  given  by  eq  (12) .  It  is  therefore  necessary, 
as  discussed  in  section  2,  to  evaluate  the  integral  of  the  force 
functions  which  is  given  below. 

In  the  present  case,  in  view  of  the  Plemelj  relation,  as 
given  by  eq  (B.42),  the  force  functions  are  directly  given  by  eqs 
(146)  and  (147)  .  The  integral  of  the  force  functions  can  now  be 
easily  carried  out.  The  result  is 

^00  0 

IF  =  J      FA(y)dy  +  FB(y)dy  =  1  ,  (156) 

^  ^00 

where  we  have  to  take  the  limit  y    =  «. 

00 

In  deriving  eq   (156) ,  we  have  accounted  for  the  fact  that  in 

A 

the  definition  of  H  (z)  in  eq  (146)  ,  the  cut  is  on  the  positive 
real  axis.  Hence,  for  the  discontinuity  in  H  (y)  at  the  real  axis, 
we  have  used  the  relations 

i  =  exp(TTL/2)    ,   -L  =  exp  (3TTL/2 )  ;  (157) 

ln(y    +  ip)  =  In  y    ,   ln(y    -  up)  =  In  y    +  2ttl;  (158) 

w  00  00  00  00 

where  p  has  a  positive  real  part  and  mod(p/yOT)   «  1. 

Similarly,  for  H  (z)  in  eq  (146)  ,  the  cut  is  on  the  negative 
side  of  the  real  axis.   In  this  case,  therefore, 

l  =  exp(TTL/2)    ,   -I  =  exp(-TTL/2)  ;  (159) 

and 
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ln(-y    +  ip)  =  In  y    +  ttl   ,   ln(-y    -  ip)  =  In  y    -  til;  (160) 

v     J  00  oo  N-*oo  J  00  x  ' 

We  find  from  eq   (158)    that  the  force  functions  given  by  eqs 
(148)   and  (149)   do  not  satisfy  the  constraint  given  by  eq  (12) .  We 
therefore  add  the  following  force  functions,   which  are  solutions 
of  the  homogeneous  Hilbert  equation,   that  is,   eqs   (143)   and  (144) 
with  their  RHS  equal  to  zero. 

F^(Y)  =  -e/7r(y2+e2)  (161) 

and 

Fg(y)  =  -e/7r(y2+e2)    ,  (162) 

where,  in  both  eqs  (161)  and  (162)  ,  the  limit  e  =  +  0  has  to  be 
taken.  The  two  force  functions  have  a  delta  function  dependence  on 
y.  In  the  limit  e  =  0,  they  approach  zero  for  all  values  of  y 
except  for  y  =  0.  Their  integral  is  finite,  and  they  make  a 
nonvanishing  contribution  to  the  H-f unctions.  The  corresponding 
H-functions  which  represent  the  homogeneous  solution  of  the 
Hilbert  equation,  are  given  below. 

Ho<z>  =    2¥T  z_1  <163> 

and 

Ho<z>  -   WZ  z_1  •  <164> 

The   integrals   of  the   force   functions  as  given  by  eqs  (161) 
and  (162)   can  be  evaluated  as  follows 

IF(J  =     Lim       [    JYo°  F^(y)   dy  +        J°  F^(y)   dy  J.  (165) 

Yoo=   00  0  ~yoo 
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Using  eqs  (157)  —  (160)  ,  and  taking  the  limit  e  =  +  0,  and  = 
oo,  we  obtain 


In  writing  the  force  functions  given  by  eqs    (161)    and  (162) 

A  B 

the  total  integral  has  been  divided  equally  between  F  and  F  but 
this  division  is  quite  arbitrary.  Since  both  the  force  functions 
are  applied  at  the  same  point,  that  is,  y  =  0,  what  matters  is 
their  total  resultant.  Thus,  in  general,  the  coefficient  in  eq 
(161)  could  be  taken  to  be  any  number — say  £  instead  of  0.5, 
provided  the  coefficient  in  eq  (162)    is  taken  to  be  l-£. 


From  eqs  (156)  and  (166)  we  see  that  if  the  force  functions 
given  by  eqs  (161)  and  (162)  are  added  to  those  given  by  eqs  (143) 
and  (144)  ,  the  constraint  given  by  eq  (12)  will  be  satisfied. 
These  forces  will  contribute  zero  stress  at  the  free  surface, 
since  they  are  solution  of  the  homogeneous  equation.  Thus  the 
final  solution  which  will  satisfy  eq  (12)  and  all  the  other 
boundary  conditions  is  given  by 

A   R  A   R  A  R 

HA'*(Z)    =     h£'*(z)    +     H^'tJ(z);  (167) 

A  B 

where  the  particular  solutions  H  '    (z)   are  given  by  eqs   (153)  and 

P  A  B 

(154)  and  the  homogeneous  solutions  HQ'  (z)  by  eqs  (163)  and 
(164) .  Although  the  two  expressions  given  by  eqs  (163)  and  (164) 
look  similar,  the  former  has  to  be  evaluated  by  taking  a  cut  on 
the  positive  half  of  the  real  axis  and  the  latter  by  taking  the 
cut  on  the  negative  half  of  the  real  axis. 

The  stress-1  Green's  functions  are  thus  given  by  eq   (144)  in 
the  UHP  and  eq   (14  5)    in  the  LHP  with  the  H-f unctions  given  by  eq 
(167).   The  displacement  Green's   functions  are  given  by  eqs  (142) 
and    (143)    in   terms    of    the   U-f unctions   which,    according   to  eqs 
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(123)  and  (124),  are  obtained  by  integrating  the  H  functions.  The 
result  is 


ljA(z)  =  "  2Wr[v(z'~Lp)  "  V<Z'LP") 

+  d  ^v(z,-Lp)   -  v(z,tp)J]   +  JL_ln(z)  (168) 


and 

.B 


uB(z)  =  "  t?T  [V(Z'L^)   "  v(z,-ip)j   +  _L_ln(z),  (169) 


where 


l   r  2 

v(z,p)  =  -  0.5   (In  z)     +  In  z  ln(z-p)   -  In  p  ln(z-p) 

"  I  ~h  (  p/z)n    1       (for  lz/pl   >  1  )  (170) 
n    n  J 

and 

v(z,p)  =  p  In  p  In  z  +  In  z  ln(p-z)   -  In  p  ln(p-z) 

+  X  -^2   (  z/p)n  (for  |z/p|   <  1  )  , 

n    n  -1 


(171) 


and  where  the  sum  is  over  all  positive  integers  n  (n— 1,2,..).  The 
last  logarithmic  term  in  eqs  (170)  and  (171)  corresponds  to  the 
homogeneous  solution. 

The    stress-2    Green's    functions    can    be    obtained  from  the 

displacement    Green's    functions    by    using    the    method  given  in 

section  2  or  directly  from  eqs  (125)  and  (126)  .  These  are  given 
below. 
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i.  x  in  UHP  : 


TA(x,x')  ^Re        -1   +    — *   * 

2    TT  L    (Z    -   p)  (Z    -   P    )  -I 

+  [  HA[(z+e)/L]  +  HA[-(z  +  c)/l] 


"  -^f-   [  HA[-z/L]   +  HA[z/L] 


+     LdA|^  HB[Z/L]    +  HB[-z/L]J  (172) 


ii.  x  in  LHP  : 


TB(x,x')   =  -  •^B-  Re   [  ^_  1  _^  J   +  LdB[^  HA[-z/l]   +  HA[z/L] 


+  — ^-JhB[-(z+c) /L]   +  HB[(z+e)/L 


+  HB[z/l]   -  HB[-z/l]  J  (173) 


The     calculated     values     of     the     displacement     field     u(x) , 

stress-1  component  t^^x)    an<3  the  stress-2   component  r32(x)  have 

been  shown  as  3-D  plots  in  figure  1-2  through  1-4  respectively.  In 

these    calculations    d,    and    d_    have    been    taken    to    0.8    and  0.2 

A  B 

respectively.  The  unit  line  force  has  been  applied  at  x'  =  (3, 
1.5).  The  vector  x'  is  chosen  so  that  it  does  not  fall  at  any 
point  of  the  mesh  where  the  Green's  function  is  calculated.  This 
is  to  avoid  the  characteristic  singularity  in  the  Green's  function 
at  x=x' .  We  see  from  figure  1-2  through  1-4  that  t31  is  0  at  x  =0 
and  u(x)  and  x22  are  continuous  at  x2=0  whereas  Z21  ^s 
discontinuous  at  *2=0. 
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An  important  point  concerning  the  numerical  calculations:  In 
deriving  the  expression  for  H  (z)  ,  we  have  taken  the  cut  on  the 
positive  half  of  the  real  axis.  This  implies  that  the  argument  of 
the  variable  changes  by  2ttl  when  it  crosses  the  real  axis  on  the 
positive  side.  This  has  to  be  properly  accounted  for  in  the 
computer  programming.  For  example  the  complex  log  should  have  the 
following  values  in  the  limit  c  =  +0  where  y  is  real  and  positive. 

ln(y  +  lg)  =  In  y  ;  ln(y  -lc)  =  In  y  +    27tl  (174) 

and 

ln(-y  +  lc)  =  In  y  +  til   ;  ln(-y  -  lc)  =  In  y  +  til.  (175) 

On  the  other  hand,  in  normal  Fortran  on  most  computers,  the  branch 
values  are  taken  as 

ln(y  +  lc)  =  In  y  ,   ln(y  -lc)  =  In  y;  (176) 

and 

ln(-y  +lc)  =  In  y  +  til  ,  ln(-y  -  lc)  =  In  y  -  til  .  (177) 

In  the  derivation  of  H  (z)  ,  the  cut  has  been  taken  on  the 
negative  side  of  the  real  axis.  In  this  case  the  argument  of  the 
variable  would  change  from  til  to  -7tl  as  it  goes  across  the  cut. 
The  branch  values  of  the  log  function  for  this  cut  will  be  the 
same  as  given  by  eqs  (176)  and  (177)  . 

In  the  numerical  calculations  for  the  H  functions,  it  is 
therefore  necessary  to  program  the  computer  so  that  the 
appropriate  branch  values  are  taken,  that  is,  as  given  by  eqs 
(174)  and  (175)  for  the  HA- function  and  by  eqs  (176)  and  (177)  for 
the  H  -function. 
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6.  Application  to  a  Cubic  Crystal  Containing  a  Z-5  Grain  Boundary 
— Plane-Strain  Problem 

As  a  further  illustration  of  the  formulation  given  in  section 
3,  we  apply  it  to  calculating  the  SI  Green's  function  for  a  cubic 
crystal  containing  a  Z-5  tilt  boundary.  We  have  chosen  this 
particular  case  for  illustration  because  of  the  current  interest 
in  the  structure  of  grain  boundaries  and  because  its  high  symmetry 
simplifies  the  calculations  considerably.  The  SI  Green's  function 
as  calculated  in  this  section  can  be  applied  to  a  variety  of 
problems  concerning  the  elastic  properties  of  the  grain  boundaries 
such  as  dislocation  pile-up  and  crack  propagation. 

A  tilt  boundary  in  a  solid  can  be  visualized  as  a  composite 
solid  in  which  two  half  solids  are  welded  together  along  the 
boundary  and  the  crystallographic  axes  of  the  solids  are  oriented 
with  respect  to  each  other.  In  particular,  the  S-5  boundary  in  a 
cubic  solid  can  be  modelled  as  follows  (see,  for  example, 
[24-25]).  Consider  two  pieces  of  the  same  cubic  material  with 
their  crystallographic  axes  parallel  to  the  coordinate  axes. 
Rotate  one  of  the  pieces  about  the  Z-axis  by  an  angle  9  = 
tan  1(3/5)  while  keeping  the  other  piece  fixed.  Weld  the  two 
pieces  together  so  that  they  form  a  planar  interface.  This 
interface  will  be  the  S-5  grain  boundary. 

The  cubic  solid,  which  we  have  chosen  for  our  calculations, 
is  stainless  steel.  The  Green's  functions  for  this  solid  without 
the  free  surface  have  been  calculated  in  I .  In  this  paper  we  use 
the  same  model  as  described  in  I.  We  label  the  solid  which  has  its 
crystallographic  axes  parallel  to  the  coordinate  axes  as  the  solid 
A,  that  is,  in  the  UHP.  The  solid,  which  has  been  rotated  as 
described  in  the  preceding  paragraph,  is  labelled  as  the  solid  B, 
that  is,  in  the  LHP. 
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The  values    [18]    of  the  elastic  constants  of  stainless  steel 

are  taken  as  follows   (in  units  of  c..):   c, ,  =  2.2,   c,  _  =  1.3,  cMM 

44'       11  12  44 

=  1.0.  These  will  be  the  elastic  constants  of  the  solid  A  in  our 
model.  The  elastic  constants  of  the  solid  B,  that  is,  the  rotated 
solid,  can  be  obtained  by  using  the  rotation  transformation  law 
for  the  fourth-rank  tensors  as  given  below: 

B  A 
cijkl  =  Sii'Sjj,Skk/Sll'  ci'j'k'l'    '  (178) 

where  S  is  the  matrix  of  rotation.  Its  elements  are 


Sll  =  S22  =  COS  9  <179> 

and 

S12  =  -  S21  =  Sin  6  ,  (180) 

where,  for  a  Z-5  boundary,  8  =  tan  1(3/5). 

The  calculation  of  various  parameters  for  the  present  model 
has  been  described  in  I.  Here  we  only  give  the  values  of  those 
parameters  which  are  needed  for  the  present  calculation.  As  given 
in  I,  in  the  present  case,  a  and  /3  which  label  the  roots  of  eq 
(A.  19),  take  only  the  values  1  and  2.  The  third  root,  a  =  3, 
corresponds  to  the  antiplane-strain  mode  which  has  already  been 
discussed  in  the  preceding  section. 

The  most  basic  parameters  required  for  the  Green's  function 
A  B 

calculation  are  pft'  ,  which  are  the  roots  of  eq  (A.  19),  and  the 
matrices  y(a)  and  g;(a)  as  defined  by  eqs  (A.  17)  and  (A.  15) 
respectively.  The  expressions  for  these  parameters  are  given 
below.  Other  matrices  can  be  easily  calculated  in  terms  of  these 
matrices  by  using  various  relations  given  in  section  2.  In  the 
expressions  given  below  the  superscripts  A  and  B  are  not  shown  for 
notational  brevity. 
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Pl  =   [(1  -  K2)/   2]1/2  +  L[(l  +  K2)/  2] 


1/2 


p2  =  -   [(1  -  K2)/   2]1/2  +   L[(l  +  K2)/   2]1/2  , 


Z      (   1  +  C  D  2   )  . 

k  (a) 

k     (a)  =  —  Z     8  d 
*21^   '             a  p0  pa  ' 

*22(a) 

za  <  C  +    Pa2  )  . 

all(a) 

Za  P«<  1  *  p0  +  t  Pa"  >  ' 

al2(a) 

Z     [  (   1  -  /3n)p  2  +  C   ]  / 

<r21(a> 

Z      (   -1  +  B„  -  C  P  2   )  , 

a22(a) 

cAA   Z     p   (  C2  ~  £  2  +  /3_ 
44     a  *av  s        1  0  '0 

2 

) 

a 


where  a  =  1  or  2 ,  and 

Z     =     -  ^  [1  -  K4]"1/2  , 
«  4cllPa 

K2  =  1  +   (5Q  /  C)  O0  +  0.5  5Q)  , 

^  =  Cll  /  C44  ' 

*0  =   <C12  +  C44^C44  ' 
and 

50  =  (C11  "  C12  "  2C44)/C44  V 

The  point  of  application  of  the  unit  force  is  taken  to  be  the 
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same  as  in  the  previous  section,  that  is,  at  x'  =  (3.0,1.5).  As 
mentioned  in  the  previous  section,  for  numerical  convenience,  the 
value  of  the  vector  x'  should  be  chosen  so  that  it  does  not  fall 
at  any  point  of  the  mesh  where  the  Green's  function  is  calculated. 
This  is  to  avoid  the  characteristic  singularity  in  the  Green's 
function  at  x  =  x' . 

The  calculated  values  of  the  the  displacement  u^x)    and  the 

stress  components  z^ix)  ,  "^l^-^  =  x\2  ^  '  and  T22^-^  nave  been 
shown  as  3D  plots  in  figures  1-5  through  1-8  respectively.  The 
plot  for  u2 (x)  is  similar  to  that  of  u1(x)  and,  therefore,  has  not 
been  given. 

We  see  from  figure  1-5,  1-7,  and  1-8  that  u1(x),  t12  (x)  and 
t22(x)  are  continuous  across  the  interface  at  x2  =  0  as  required 
by  the  boundary  conditions  given  by  eqs  (2)  and  (3)  .  The  stress 
component     x\\^^    i-s  n°t  continuous  across  the  interface  as  shown 

in  figure  1-6.  The  components  T21^-^  anc^  Tll^-^  are  zero  a^  ^e 
free  surface,  that  is,  at  x1  =  0  as  required  by  the  boundary 
condition  given  by  eq  (4)   and  (5) . 

The  singularities  in  the  stresses  will  be  discussed  in  the 
following  paper  on  the  generalized  plane-strain  problem.  As 
mentioned  in  section  3D,  the  singularities  in  H(z)  arise  from  the 
homogeneous  part  of  eq  (36)  and,  in  general,  do  not  depend  upon 
the  RHS  of  eq  (36)  .  As  we  shall  see  in  the  following  paper,  the 
matrix  M(q)  is  exactly  the  same  in  the  present  plane-strain  case 
and  the  generalized  plane-strain  case. 

It  may  be  emphasized  again  that,  in  using  the  formulation  of 
this  paper  for  numerical  calculations,  the  values  of  the  various 
functions  such  as  ln(z)  and  zL(^  over  the  upper  and  lower  branches 
of  the  cut  have  to  evaluated  carefully  according  to  the  scheme 
given  by  eqs  (174)-(177). 
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7.  Summary  of  Main  Results 

In  this  paper  we  have  obtained  the  displacement  Green's 
functions  and  the  stress  Green's  functions  for  a  general 
anisotropic  bimaterial  composite  solid  which  has  a  plane  interface 
and  a  free  surface  normal  to  the  interface.  These  are  applicable 
to  a  variety  of  plane-strain  problems  and  also  account  for 
displacements  and  forces  normal  to  the  plane. 

We  find  that  the  stress  field  is  singular  near  the 
intersection  of  the  free  surface  and  the  interface.  The 
singularities  are  of  the  type  r  11  and  [ln(r)]n  where  7)  is  a 
positive  number  between  0  and  1  and  n  is  a  positive  integer 
between  0  and  6.  The  actual  values  of  7)  and  n  and  the  coefficient, 
i.e.,  the  weight  of  the  singularities,  depend  upon  the  material 
properties  of  the  solids.  In  addition  to  the  singularity,  an 
oscillatory  behavior  of  cos  or  sin[ln(r)]  type,  may  also  be 
present  in  certain  cases. 

A  method  of  solution  of  the  generalized  inhomogeneous  vector 
Hilbert  problem  has  been  given.  The  method  is  simple  and 
convenient  for  applications  to  different  problems  in  the  stress 
analysis  of  materials  involving  cracks  and  interfaces  in  which  the 
inhomogeneity  and  the  nonsingular  part  of  the  Hilbert' s  kernel  are 
of  the  form  l/(y  -  p)  where  y  is  a  real  variable  and  p  is  complex. 
For  this  purpose,  we  have  introduced  a  complex  transform  which  is 
of  the  form  yL<^  0,5  where  l  =  v'(-l)  and  q  is  a  real  variable.  This 
transform  seems  to  be  particularly  interesting,  because  it  is  an 
eigenf unction  of  the  Hilbert' s  kernel  and  is  orthogonal  over  both 
variables  y  and  q  in  the  domain  0  ^  y  <  «  and  -«  <  q  <  w. 
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Appendix  I -A 


Elastic  Green's  Function  for  an  Infinite  Bimaterial  Composite 

In  this  appendix  we  quote  the  expressions  for  the  elastic 
Green's  function  for  a  bimaterial  composite  as  obtained  in  I  and 
used  for  calculation  of  Green's  function  for  a  cracked  solid  in 
II.  The  expressions  given  below  have  been  taken  from  II  and  are  in 
a  form  which  is  more  convenient  for  present  application.  The 
coordinate  system  used  for  these  expressions  is  the  same  as  shown 
in  figure  1-1.  The  indices  a  and  /3  take  the  values  1,   2,  or  3. 

The  Green's  function  given  below  is  the  displacement  Green's 
function  G(x,x'),  which  gives  the  displacement  field  at  the  point 
x  when  a  unit  line  force  is  applied  at  the  point  x' .  As  in  I  and 
II,  we  distinguish  between  four  different  cases  corresponding  to  x 
and  x'  being  in  UHP  and  LHP. 

1.  x  and  x'  both  in  UHP  (  x    *  0  ;  x  '  *  0  )  : 


a 


(A.l) 


a/3 


2  . 


x  in  LHP, 


x'  in  UHP  (  x2 


s  0 


•  x  ' 
'  X2 


2:  0   )  : 


Gq   (x,x')  = 


1  ~B,   v    ^11  -,    ,-B     -A  , 

^  2_    *  (0°  2f?    ln(zoT  P/3  > 


(A. 2) 


a/3 
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3.  x  in  UHP,  x'   in  LHP  (  x2  £  0  ;  x2'  ^  0  )  : 

GqB(x,x')   =  .  1 rA(cx)   q"1  ln(z^-  pA  )    ,  (A. 3) 

a/3 

4 .  x  and  x7  both  in  LHP  (  x2  ^  0  ;  x2 '  £  0  ) : 

GBB(x,x')  =  -  I  £  JB(a)   In (5*-  pB) 

a 

1  V     "B,   v    ^IV  ,    ,-B  B. 
+  n  )      I   (a)   2|3     ln(za~  <V  '  (A-4) 

a/3 

where 

A,  B                 A.B  /T.  ... 

V      =  Xl  +  V      X2  '                                                                                    (A»  5) 

pA,B  =  xl   +  pA/B  xl   ,  (A. 6) 

Qa  =  M  [  aA(oc)   -  eB(iB)_1  rA(a)]    ,  (A. 7) 

Q*1  =  N  [  aA(a)   ~  /(is1"1  IA(«)]    »  (A.  8 ) 

Ql11  =  M  [  ciB(cO   -       (Ig)"1  *V)]    ,  (A.  9) 

91V  =  n  [  eB(«)  -  ^s^s1"1  *B<a>]  '  (A-10> 


,  A.  -1    T  ~B    ,-B.  -1       A,  A. -ll"1  „„, 
M  =   (Is)        [  2S   ^s*     "  2s(2S)     J  '  (A.  11) 


and 
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,-B.  -1    T  -B    ,-B.  -1       A,  A.-ll 


-1 


(A. 12) 


(Note.  The  following  formulas  are  valid  for  either  superscript,  A 
or  B,  corresponding  to  the  solid  in  UHP  or  LHP,  respectively.) 


Is  =  I  ' 
a 


(A.  13) 


Zs  =  I  2t«)  / 


(A. 14) 


a(a)  =  L1Z(a)  r(oc)  , 


(A. 15) 


II 


Lij«*>   =  Ci2kl  +  Pa  Ci2k2  ' 


(A.  16) 


K  •  • (a) 


rii<q2> 


aq. 


<p«  -  p«'  n 

a*/3 


(for  q2  =  q    pj  , 


(A. 17) 


• (g)   is  the  ij  cofactor  of  the  Christoffel  matrix  A  defined  by 


Aij    (3)   =  cikjl  qk  qx 


(A. 18) 


qn  and  q_  are  components  of  the  wave  vector  g  and  p  is  obtained 
such  that  q2  =  pa  q1  is  a  root  of  the  equation, 


I |A(q) ||=0, 


(A. 19) 


|A|  denotes  the  determinant  of  the  matrix  A  and  a  is  the 
coefficient  of  q.     in  | |A| | ,  and 


Im  (p  )   >  0  . 


(A. 20) 
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In    analogy    with    eqs     (A. 13)     and     (A. 14),     we    define  the 
following  quantities  as  sums  over  a: 


<4 


,11, in, iv  =  y  Qi, ii, in, iv  _  (A  21) 


a 


Several  useful  relations  between  the  various  parameters  of 
the  Green's  function  (such  as  j,  cr,  Q  )  have  been  given  in 
Appendix  A  of  II.  Two  more  relations  which  involve  the  constants 
introduced  in  this  paper  are 

AD  r  =  0  ,  (A. 22) 


where 


n(a)/Pa  =  -v(cc)  =  J(a)    ,  (A. 23) 


*ik  =  Cilkl  +  PcJ  Cilk2  +  Ci2kl>    +  Pa  Ci2k2    '  (A. 24 ) 
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Appendix  I-B 


Solution  of  Generalized  Inhomogeneous  Vector  Hilbert  Problem 


In  this  appendix  we  give  a  solution  of  the  generalized 
inhomogeneous  vector  Hilbert  problem.  The  corresponding  scalar 
problem,  its  properties  and  solution  have  been  discussed  in  detail 
in  the  excellent  treatise  by  Muskhelishvili  [20] .  The  vector 
problem  has  been  discussed  by  Vekua  [21]  .  Our  method  of  solution 
as  given  below  is  different  from  those  given  in  [20,21]  and  is 
based  upon  the  use  of  a  complex  transform.  This  method  should  be 
particularly  convenient  for  those  types  of  kernels  and 
inhomogeneities  which  occur  in  problems  of  stress  analysis  in 
solids  containing  surfaces,   interfaces,  and  cracks. 

First  we  consider  the  scalar  problem  as  defined  below: 

D  H(y+LC)   +  6  H(y-LC)   +  £  Ca  H(aay)   =  f(y),  (B.l) 

a 

where  D,  Ca  and  a^  are  complex  numbers,  a  is  a  summation  index 
which  takes  an  arbitrary  but  finite  range  of  values,  t  and  y  are 
real  variables  which  vary  between  0  and  »,  c  is  a  small  positive 
number  which  approaches  zero  in  the  limit,  f  (y)  is  a  known 
function,  and  H(z)   for  any  complex  z,   is  defined  as 


CO 


H(z)   =  1 


2ttl 
0 


F(t)     dt  .  (B.2) 


t  -  z 


The  integral  in  eg  (B.2)  has  to  be  interpreted  as  the 
principal  value  integral  for  z  on  the  real  axis  when  the  integral 
is   singular.    The   first   two   terms    in   eq    (B.l)    correspond   to  the 
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principal  value  of  the  integral  in  the  definition  of  the  H- 
f unction.  The  function  F(t)  in  eq  (B.2)  is  the  unknown  function 
which  has  to  be  determined  from  eq  (B.l). 


In  eq  (B.l),  f (y)  is  the  inhomogeneity ,  which  is  zero  for  the 
homogeneous  Hilbert  problem.  If  C  =  0,  this  equation  will  have  no 
nonsingular  term  in  the  kernel  and  is  then  called  the  special  (not 
generalized)  Hilbert  problem  (the  prefix  'special'  is  usually 
omitted) .  The  particular  form  of  the  nonsingular  kernel  as  given 
in  eq  (B.l)  is  of  the  type  1/ (t-ay)  where  a  is  complex,  and  is  of 
special  interest  in  problems  concerning  stress  analysis  in  solids 
containing  surfaces,   interfaces,  and  cracks. 

We  solve  eq  (B.l)  by  using  a  complex  transform  method.  For 
this  purpose  we  introduce  the  following  function  of  two  real 
variables,  q  and  y  (0  ^  y  ^  »  and  -co  <  q  s  co)  : 

T(q,y)   =  yLq-°-5   .  (B.3) 


By  using  the  values  of  the  branch  points  as  given  in  Appendix  I-C, 
we  can  verify  that  T(q,y)  is  an  eigenf unction  of  the  homogeneous 
Hilbert  problem  given  by  the  LHS  of  eq  (B.l).  Further,  the 
functions  T(q,y)  are  orthogonal  in  q  as  well  as  y  space  as  given 
below: 


.00 


T(q,y)   T(q',y)   dy  =  2tt  5  (q-q'  ) 


(B.4) 


0J 


and 


.00 


T(q,y)   T(q  ,y')   dq  =  2tt  5  (y-y' ) 


(B.5) 


-co 
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These   orthogonality  relations   can   be  verified  by  using  the 
following  standard  definition  of  the  delta  function: 


5(y)  = 


2TT 


exp(iqy)  dq 


(B.6) 


The  integral  in  eq  (B.4)  can  be  evaluated  by  using  the 
substitution  x  =  In  y.  In  writing  eq  (B.5),  we  have  used  the 
following  property  of  the  delta  function  of  an  arbitrary 
dif f erentiable  function  p(y)   of  y: 


5[P(Y)]  =  S(y-y0)    [  Hv^l   1       '  (B*7) 

where  y  =  yQ  is  a  real  root  of  the  equation  p(y)  =  0. 

We  now  write  F(t)   in  terms  of  its  transform  v(q)  as 


n  1  _  _ 

F(t)    =  v(q)    tiq"         dq   .  (B.8) 

—  00 

Substituting  for  F(t)   from  eq   (B.8)    in  eq   (B.2),  and  carrying  out 

the  integration  over  y  as  shown  in  Appendix  I-C  (see  eq  C.7),  we 
obtain  the  following  result  for  H(z) : 

H(Z)    =        f     V(q)    ziq~°'5  dq   ,  (B.9) 
-coJ 

where 

V(q)   =  v(q)/E(q)  (B.10) 

and 
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E(q)   =  l/[l+exp(-27rq)  ] 


(B. 11) 


As  described  in  Appendix  I-C,  in  deriving  eq  (C.7),  we  have 
taken  a  cut  on  the  positive  real  axis  (see  figure  1-9) .  The  values 
of  argument  just  above  and  just  below  the  real  axis  are  therefore 
taken  to  be  0  and  2tt,  respectively,  as  given  in  eqs  (C.3)  and 
(C.4).  Thus,  we  get  the  following,   in  the  limit  e  =  +0, 

00  , 

H(y+ie)   =  V(q)   ylc5"0*5  dq  (B.12) 

-or 

and 

00  . 

H(y-te)   =     -[     V(q)  exp(-27rq)   ylc3"0*5  dq  .  (B.13) 


In  order  to  solve  eq  (B.l),  we  use  the  representation  of  H(z) 
as  given  by  eq  (B.9)  and  use  eqs  (B.12)  and  (B.13)  for  the  first 
two  terms.  Then  we  multiply  both  sides  of  eq  (B.l)  by  T(q,y)  , 
i.e.,  y  LC*  ®'  ,  and  integrate  over  y  from  0  to  oo.  Using  the 
orthogonality  relation  as  given  by  eq  (B.4),  we  obtain  the  result 

V(q)   =   [M(q) j"1  N(q)  (B.14) 


where 


M(q) 


=   [  D  -  D  exp(-27rq)   +  £  Cft  a^q  °*5  ] 


(B. 15) 


and 


00 


x=  t   \      -Lq-0.5  , 
f(y)   y  dy 


(B. 16) 


0J 
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In  many  cases  of  practical  interest  such  as  stress  analysis 
in  solids,  the  function  f (y)   is  of  the  form 


f (Y)  =  X  R«  /  [y  "  Pa]    '  (B'17) 


a 

where  pa  is  complex.  In  such  cases,  the  integral  on  the  RHS  of  eq 
(B.16)  can  be  easily  obtained  by  using  the  same  contour  as  used 
for  deriving  eq  (C.7)  in  Appendix  I-C.  The  result,  which  can  be 
obtained  by  replacing  q  by  -q  (not  i  by  -l)  in  eq  (C.7)  is  given 
below: 


=  Lexp(-27iq)   j  -iq-0.5 

E(q)       L    a   Vla'  v  ' 


Equation   (B.14)    along  with  eqs   (B.8),    (B.9)    and   (B.10)  gives 

the    particular    solution    of    eq    (B.l).    In    order    to    obtain  the 

general    solution,    we   have   to   add   a   solution   of   the  homogeneous 

part    of    eq    (B.l)     to    the    particular    solution.    The  homogeneous 

solution  can  be  written  in  terms  of  Q  as 

r 


lQ  -0.5 

H  (z)   =  £  h(Q  )    z     r  ,  (B.19) 

r 


where  h(Qr)  are  arbitrary  constants;  and  Q  ,  which  may  be  complex, 
is  a  solution  of  the  following  equation: 


M(Qr)   =  0   .  (B.20) 

The  arbitrary  constants  can  be  determined  with  the  help  of  the 
boundary  conditions  on  F(t). 

Now  we  consider  the  vector  problem  which  arises  if  there  is 
more    than    one    unknown    function    F(t)    and,    of    course,    an  equal 
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number  of  equations.  In  this  case  we  regard  eq  (B.l)  as  a  matrix 
equation.  If  there  are  n  unknown  functions  and  n  equations,  then  D 
and  in  eq  (B.l)  will  be  n  x  n  square  matrices,  whereas  H(z)  , 
f  (y)  will  be  n  x  1  column  matrices  or  n-component  vectors. 

Similarly,  M(q)  and  N(q)  will  be  n  x  n  and  n  x  1  matrices 
respectively.  Equation  (B.20)  will  become  a  determinantal  equation 
and  h(Qr)  will  be  an  eigenvector  of  M(Qr)  corresponding  to  a  zero 
eigenvalue.  One  of  the  components  of  this  eigenvector  will  remain 
arbitrary,  which  has  to  be  determined  from  the  boundary 
conditions.  With  these  modifications,  the  solution  given  above 
remains  valid  for  the  vector  Hilbert  problem. 

As  a  special  case  of  the  general  method,   we  solve  the  vector 

Hilbert  problem  for  two  unknown  functions.    This  solution  is  used 

in  this  paper  (section  2) ,  where  the  unknown  functions  are  the  two 

A  A 

force  functions  F  and  F  .  As  in  section  2 ,  we  label  the  unknown 
functions  with  A  and  B  and  also  label  the  corresponding  matrix 
elements  of  eq  (B.l)  by  the  same  indices.  Equation  (B.l)  now 
becomes  a  set  of  two  simultaneous  equations  as  given  below: 


D 


AA 


AA 


a 


(B.21) 


a 


and 


+  D 


BB 


BB 


H  (y-LC) 


(B. 22) 


a 
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The  H-f unctions  are  now  defined  as 


"a  FA(t) 


t  -  z 


dt 


(B.23) 


and 


HB(Z)   =  1 


2TII 


-oo" 


°  FB(t) 


t  -  z 


dt 


(B.24) 


As    in    eq    (B.8),    the    force    functions    are    written    in  the 
following  form: 


FA(t)    =       f  vA(q)   t^"0-5  dq 


and 


—  00 


FB(t)    =        f  vB(q)  |t|i(J-0-5 

—  00 


dq 


(B.25) 


(B.26) 


The  corresponding  expressions  for  the  H-functions  are 


HA(Z)  -- 


=       J     Vft(q)  z 


iq-0 . 5 


dq 


(B.27) 


and 


HB(Z)   =       f  VB(q)   z^5  dq  . 

—  00 


(B.28) 


Using  eq    (B.25),    the   integration  for  H  (z)    in  eq    (B.23)  is 

carried  out  as  in  eq  (C.7).  The  result  is  identical  to  that  given 

by  eqs    (B.9)    and    (B.10)    except  that  the  superscript  A  has  to  be 

added   to   the   relevant   matrices.    The   corresponding    integral  for 
g 

H  (z)  has  also  been  evaluated  in  Appendix  I-C  [see  eq  (C.13)]. 
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In  case  of  H  (z)  ,  we  have  taken  the  cut  on  the  negative  real 
axis.  The  argument  of  z  is  taken  to  be  n  just  above  the  real  axis 
and  -7i  just  below  it.  The  corresponding  values  of  the  H-f unctions, 
as  obtained  by  using  eqs  (C.ll)  and  (C.12)  are  given  below:  (in 
the  limit  c  =  +0;  -«  <  y  <  0)  : 

00  ■ 

HB(y+LC)   =     -lJ     VB(q)exp(-7rq)    |y|lc*-0-5  dq  (B.29) 


—  00 


and 


HB(y-ie)   =       l\    VB(q)exp(7rq)    |y  |         . 5  dg  #  (B.30) 

—  00 


To  solve  eqs  (B.21)  and  (B.22),  we  proceed  as  before,  except 
that  we  multiply  eq  (B.21)  by  y~L(*~0'5  and  eq  (B.22)  by  |y|"Lq"0,5 
and  then  integrate  over  y  between  the  limits  given  in  eqs  (B.21) 
and  (B.22)  to  ensure  the  orthogonality  of  the  transform.  Thus,  we 
obtain  the  following  particular  solution  for  the  force  functions 
and  the  corresponding  H-f unctions: 

VA(q)   =  VA(q)/E(q)   =  NA  +  ^   ,  (B.31) 


VB(q)   =  ivB(q)exp(-7rq)/E(q)   =  M*A  NA  +  Mj^  NB  ,  (B.32) 

where  M1  is  the  inverse  of  M.  The  matrix  elements  of  M  and  N  are 

MAA  "  DAA  "  5AA  ««P(-2"q)    +  £  ^  (aj           '  5   ,  (B.33) 

a 


M 

AB 

a 


■     I  WViq"°-5   <  (B.  34) 
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m     =  V  c     r-c  )i(J-°-5 

BA      L    aBA  *  ca; 
a 


(B.35) 


MBB  ="L[DBB  exP(-"q) 


-  DBB  exp(Trq)]   +  J  CaBB         )  iq"° '  5   f  (B.36) 
J  a 


.00 


NA(q) 


1 

27T 


_   ,   .      —  iq—  0  .  5  , 

fA(y)  y   H  dy 


(0  *  y  *  oo)    ,  (B.37) 


and 


NB(q) 


_1 

27T 


oo 


fB(y) 


y|    L(3"0-5  dy     (-oo  <  y  <  0) 


(B.38) 


Equations  (B.31)  and  (B.32),  when  substituted  into  eqs  (B.27) 
and  (B.28),  give  the  particular  solution  for  the  H-functions.  The 
homogeneous  solution  is  given  by  eq  (B.19),  where  h(Qr)  is  now 
identified  as  an  eigenvector  of  the  matrix  M(Qr)  for  a  zero 
eigenvalue  with  Qr  as  a  solution  of  the  following  determinantal 
equation, 


| |M(Qr) ||=0. 


(B. 39) 


The  final  solution  is  the  sum  of  the  particular  solution  and 
the  homogeneous  solution  as  given  below: 


rA(z)  =  r  vA(q)  z±q  0,5  dq  + 1  w z  r 

-ooJ  r 


iQ  -0.5 


(B.40) 


and 


r  r°°  -irr-n  r  ^  lQ  -0.5 

HB(z)   =  V   (q)    ziq  U*D  dq  +  £  hg(Qr)    z     r  .  (B.41) 

-ooJ  r 
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One  of  the  two  components  of  the  eigenvector  h,  (Q  )  and  h_,(Q_  )  is 
arbitrary  and  has  to  be  determined  from  the  boundary  conditions  or 
the  constraints  imposed  on  the  force  functions. 

The     force     function     can     be     obtained     directly     from  the 
H-f unctions  by  using  Plemelj  relation  [20] 

F(y)   =  H(y+te)    -  H(y-Le)    .  (B.42) 

In  the  present  case,  eq  (B.42)  may  be  easily  verified  by  using  the 
transforms  of  F  and  H  as  given  by  eqs   (B.31)   and  (B.32). 

In  many  problems   concerning   stress   analysis  in   solids,  the 

inhomogeneity    f(y)    is    of    the    form   given    by    eq  (B.17).    In  the 

present  case,  we  write  the  two  components  of  the  vector  f(y)  in 
the  form, 


fA(y)   =  I  Ra  /  [y  "  pa]  (0  "  y  "  °°)  (B*43) 


a 


and 


fB(y)  =  £  R*  /  [y  -  pjj|]  (-00  <  y  <  0)    .  (B.44) 

a 


The  values  of  the  integrals  for  NA (q)  and  N„(q)  in  eqs  (B.37) 
and  (B.38)  can  be  obtained  from  eqs  (C.7)  and  (C.13)  by  replacing 
q  by  -q  and  not  l  by  -l.  The  result  for  NA(q)  is  identical  to  that 
given  by  eq  (B.18)  except  that  the  superscript  A  has  to  be  added 
to  R    and  p     on  the  RHS.  The  result  for  ND(q)  is 


V«>   --SSE^XR2   (Pa>-tq-°-5   •  (B.45, 


71 


Appendix  I-C 


Evaluation  of  Integrals 

In  this  appendix,  we  shall  evaluate  certain  integrals  by 
using  contour  integration.  These  integrals  are  required  for 
deriving  various  formulas  in  the  paper.  We  shall  also  obtain  the 
force  function  and  its  integral  to  ensure  that  the  constraint 
given  by  eg  (12)   is  satisfied. 


required  for  calculating  the  H     function  from  the  force  function 


First,    we    shall    evaluate    the    following    integral,    which  is 
.red  for  calculate 
in  eq  (B.  2)   or  (B.  23)  . 

00  .  _ 

tiq-0 . 5 


IA(q,P)  = 


t  -  p 


dt  ,  (C.l) 


where  p  is  complex  and  q  is  real. 

To  evaluate  I   (q,p) ,    consider  the   integral  of  the  following 
complex  function  of  variable  z  : 

iq-0.5 

Z  =    —   ,  (C.2) 

z  -  p 

where  z  is  complex.  We  integrate  Z  over  the  contour  shown  in 
figure  1-9.  The  contour  consists  of  an  outer  circle  C  with  a  cut 
on  the  positive  real  axis.  The  radius  of  the  circle  approaches 
infinity  in  the  limit.  The  argument  of  z  is  taken  to  be  0  at  the 
upper  branch  (just  above  the  cut)  and  2n  at  the  lower  branch  (just 
below  the  cut) ;  that  is, 

t    =  t  +  ie  =  t  (C.3) 
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and 

t_  =  t  -  lc  =  t  exp(27TL)    ,  (C.4) 
where  t  is  the  real  part  of  z. 

Thus  we  obtain  from  Cauchy's  theorem, 


Ziq-°-5     i>  .  .A 


C     Z  -  p 


dz  +  E(q)    IA(q,p)   =  2ttl         Residues  (C.5) 


where 


E(q)   =  1  +  exp(-27rq)    .  (C.6) 

The  first  term  in  eq  (C.5),  denoting  the  contribution  of  the 
outer  circle,  vanishes  in  the  limit  when  C  extends  to  infinity. 
Since  the  integrand  has  only  one  pole  inside  the  contour  at  z  =  p, 
we  obtain  from  eq  (C.5) 


i.q-0.5 

TA,  ,  2ni  p  ^  -7\ 
1   (q'p)   =   Efq)   *  (C'7) 


In  order  to  calculate  the  value  of  H  (y+te)  and  H  (y-te) ,  we 
need  the  value  of  I(q,p)  just  above  and  just  below  the  real  axis. 
Using  eq  (C.3)  and  (C.4),  we  obtain  the  following  values  of  the 
integral  over  the  upper  and  the  lower  branches  of  the  cut: 


I   (q,y+Le)  = 


2ttl  y 


Lq-0 . 5 


E(q) 


(C.8) 


and 


I   (q,y-te)  = 


2?Tiexp (-27rq)  y 
E(q) 


(.q-0 .  5 


(C.9) 
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We  now  evaluate  the  following  integral,  which  is  required  for 
g 

calculating  H  (z)   in  eq  (B.24): 


IB(q,P)  = 


-co 


°  |+-lic3-0-5 

-L^J   dt  .  (C.10) 

t  -  p 


g 

To  evaluate  I  (q,p) ,  we  consider  the  same  function  as  defined  by 
eq  (C.2)  and  evaluate  its  integral  on  the  contour  given  in  figure 
1-10,  which  has  a  cut  on  the  negative  real  axis.  The  argument  of  z 
is  taken  to  be  n  on  the  upper  branch  and  -n  on  the  lower  branch 
so  for  -co  <  t  ^  0, 

t    =  t  +  lc  =  |t|exp(7Ti)  (C.ll) 

and 

t_  =  t  —   LC   =    1 1 1  exp  (-7TL )    .  (C.12) 


Proceeding  as  before,  we  obtain 


IB(q,p)   =  -  2"eXP(l(Vtq~0"5      •  (C.13) 


g 

The  values  of  I  (q,p)  just  above  and  below  the  real  axis  as 
calculated  by  using  eqs   (C.ll)-  (C.13)  are 


IB(q,y+tc)    =     2TIcexp(-21Tq)|y|'-q-0-5  (c  14) 


and 


i    i Lq-0 . 5 

TB.  ,  27TL    y      ^  ,~  i<=\ 

I   (q,y-LC)   =   .  (C.15) 
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We  now  evaluate  the  integral  over  q  required  for  the 
calculation  of  the  H  function.  As  we  see  from  eq  (54) ,  we  need  to 
evaluate  an  integral  of  the  form, 

CO 

P(q)   dq     zLq~0'5       ,  (C.16) 
E(q)D(q) 

—  00 

where  D(q)  is  the  determinant  of  the  matrix  M(q)  and  P(q)  is  an 
analytical  function  of  q  with  no  poles. 

To  evaluate  the  integral  in  eq  (C.16),  we  choose  a 
semicircular  contour  either  in  the  UHP  or  in  the  LHP  as  shown  in 
figure  1-11.  In  the  present  case,  the  contribution  of  the  integral 
over  the  large  semicircle  vanishes  in  the  UHP  if  mod(z)  >  1  and  in 
the  LHP  if  mod(z)  <  1.  Accordingly,  we  choose  the  contour  in  the 
UHP  if  mod(z)   >  1  and  in  the  LHP  if  mod(z)   <  1. 

Let     qn     and     Qr     denote     the     zeroes     of     E(q)      and  D(q) 

respectively.    The   integrand   in  eq    (C.16)    has  two   sets  of  poles; 

those  at  q    will  be  called  Set  I  and  those  at  Q    Set  II.  Thus 

r 

qn  =   (n  +  0.5)i  ,  (C.17) 

where  n  is  an  integer-  which  can  be  positive,  negative  or  zero  and 
are  the  roots  of  the  determinantal  equation, 

D(Qr)   =  0  .  (C.18) 

We  assume  that  all  Qr  are  distinct,  which  should  be  generally 
true.  When  they  are  not,  it  would  be  a  higher  order  pole  which  can 
be  treated  by  the  standard  methods  of  contour  integration. 
Moreover,  in  general,  some  q^  will  be  equal  to  some  Q^..  These 
second-order    poles    can    also    be    treated    by    standard  methods. 
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However,  presently  our  object  is  to  write  a  general  expression  for 
the  integral  with  the  view  to  identify  its  singularities  and  to 
obtain  an  expression  for  the  force  function  and  its  integral.  For 
this  purpose,  it  would  be  convenient  to  introduce  a  small  positive 
imaginary  part  lc  in  q  in  the  denominator  of  the  RHS  of  eq  (C.16) 
as  follows  and,   in  the  end,  take  the  limit  c  =  0. 


00 


IH(z)  = 


P(q)dq  ziq  0,5   .  (C.19) 

E(q-ie) D(q+ie) 


—  00 


The  two  sets  of  poles  in  the  integrand  in  eq   (C.19)    can  now 
be  written  as  follows: 

Set  I.   Zeroes  of  E(q) . 


q  =  qn  +  ic   .  (C.20) 


Set  II.   Zeroes  of  D(q) . 


q  =  Qr-ic   .  (C.21) 


The  allowed  values  of  qn  and  Qr  depend  upon  whether  the  contour  of 
integration  (see  figure  1-11)  has  been  chosen  to  be  in  the  UHP  or 
LHP  and  are  given  below: 

I.  Mod(z)   >  1   (Contour  in  UHP). 

n  =  0 ,  1 ,  2 ,  .  .  . 

and 

Im(Qr)   >  0.  (C.22) 

II.  Mod(z)   <  1   (Contour  in  LHP). 
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and 


Im(Qr)   <  0 


(C.23) 


The  residue  of  E(q-ie)  at  qn  is  2tt  and  that  of  D(q+ie)  at  Qr 
is  D'  (Q  ) ,  where  the  prime  denotes  differentiation  with  respect  to 
the  argument.  Thus,  we  obtain  the  following  value  of  the  integral 
by  using  Cauchy's  integration  formula: 


IH(z) 


^  D(q 


P(q  +ie)  z-(n+1>-e 


_ (q  +  2ic) 
n         v^n  ' 


iQ  -0.5+c 
P(Qr-ie  )z  r 


*  2"LY.  E(Qr-2;g)D'(Qr)  *  R<Z'C>    <  (C«24) 

r 


where  the  sums  in  eq  (C.24)  are  over  those  values  of  qn  and  Qr 
which  are  not  degenerate,  that  is,  which  represent  only  a 
single-order  pole;  and  R(z,e)  denotes  the  contribution  of  those 
values  of  n  and  r  which  represent  a  pole  of  second  or  higher  order 
in  the  integrand.  The  plus  signs  in  eq  (C.24)  are  to  be  taken  when 
mod(z)  >  1  and  negative  when  mod(z)  <  1. 

In  the  cases  considered  in  section  3A  and  3B,  D(q)   is  zero  at 
the  following  two  values  of  q  which  are  part  of  Set  I: 

qo  =  Qo  =  L/2  (C.25) 

and 

q_1  =  Q_±  =  -l/2   ,  (C.26) 

which  correspond  respectively,  to  n  =  0  and  -1  in  eq  (C.17).  These 
values  therefore  represent  second-order  poles  which  will 
contribute,    respectively,    when  the  contour   is   chosen   in  UHP  and 
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LHP.  From  eqs  (C.22)  and  (C.23)  the  pole  at  qQ  will  contribute 
when  mod(z)   >  1  and  the  one  at  will  contribute  when  mod(z)  < 

1.  We  shall  now  evaluate  R^(z,c)  at  these  poles. 

First  we  consider  the  former  case,    when  mod(z)    >   1.    We  use 
the  expansion  formulas, 


E(qQ-2L£)   =-47TLe(l  +  2nic) 


(C.27) 


and 


D(qo+  2  lc)   =  2  lc  D'(qQ)    [1  +  LCD"  (qQ)  /  D'  (qQ)  ] 


(C.28) 


Thus,  we  obtain  the  following  for  mod(z)   >  1: 

-1 


Lim  R  (z,e)  = 
e->o     u  D'  (qQ)  |_ 


P'  (q0)   +  LP(qQ)  Ln(z) 


^  P(q0)D"  (qQ)/D' (qQ)  +TrP(qo) 


(C.29) 


The  corresponding  expression  for  R(z,e)   when  mod   (z)    <  1  is 
formally  similar  to  eq  (C.29)   and  is  written  as 

Lim  R,  (z , e)  =  — - — 
e->o     1  D'  (qQ) 


P'(q_1)   +  LP(q_1)  Ln(z) 


1 

2 


P(q_1)D"  (q^/D'  (q_1)  +7TP(q_1) 


(C. 30) 


We  now  calculate  the  integral  required  for  the  force 
functions.  An  explicit  knowledge  of  these  functions  is  not 
required   for   calculating  the  Green's   functions   for  which,    as  we 
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see  from  eqs  (116)  and  (117)  ,  a  knowledge  of  H-f unctions  is 
sufficient.  However,  we  do  need  to  calculate  the  integral  of  the 
force  functions  in  order  to  satisfy  the  constraint  given  by 
equation  (12) . 

Using  the  Plemelj  relation,  as  given  by  eq  (B.42),  we  see 
that  the  two  force  functions  can  be  written  in  terms  of  an 
integral  of  the  form, 


CO 


E(q)  P(q)dq  | y | iq  °*5  .  (C.31) 
E(q-ie)D(q+ie) 


—  00 


In  case  of  the  force  function  in  the  LHP,  the  numerator  of 
the  integrand  in  eq  (C.31)  will  have  a  factor  exp(Trq)  which  can  be 
easily  included  in  the  result  derived  below.  The  integral  in  eq 
(C.31)  can  be  evaluated  by  using  the  same  method  as  described 
earlier.  The  result  is 


P(qn+ie)|yf  (n+1)"C 


I    (y)=  -2TTG 

D(q  +  2ic) 


+  2711 


E(Qr-Lc)P(Qr-ic  ) |y| 
E(Qr-2(.e)D'  (Qr) 


iQ  -0.5+e 
r 


+  R(y/p,e) ,  (C.32) 


where  the  allowed  values  of  n  and  r  are  defined  by  eqs  (C.22)  and 
(C.23),  and  R(y/p,e)  accounts  for  any  second  order  poles  as  the 
corresponding  R  term  in  eq  (C.24).  The  complex  variable  p  is 
contained  in  the  function  P(q) . 

The  first  term  on  the  RHS  of  eq    (C.32)    (along  with  the  last 
term  where   it   exists)    gives  the  particular   solution  of   eqs  (25) 
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and  (26)  .  The  second  term  is  a  solution  of  the  homogeneous 
equation  since  D(Qr)  is  zero  for  all  Qr»  The  first  term  on  the  RHS 
of  eq  (C.32)  is  vanishingly  small  in  the  limit  e  =  0.  However, 
this  term  can-not  be  neglected  because  it  gives  a  finite 
contribution  to  H(z)  ,  as  we  can  see  from  eqs  (C.19)  and  (C.24). 
The  1/y  term  which  corresponds  to  n  =  0  or  q  =  i/2 ,  is  actually  a 
delta-function  term  which  is  zero  for  all  values  of  y  except  y=0. 
This  is  similar  to  the  term  in  eqs  (158)  and  (159)  for  the 
antiplane-strain  mode. 

The  integral  of  the  force  function  in  the  q  representation 
can  be  directly  obtained  from  eq  (C.32).  We  write  the  total 
integral  of  the  force  function  as  a  sum  of  those  for  the  UHP  and 
the  LHP  as 


it  =  - 


+  I 


B 


(C. 33) 


where 


ip(y)  <Jy 


oo 


—  00 


PA(q)dq 


D(q) (iq+0.5) 


iq+0 . 5 


00 


iq+0. 5 


(C. 34) 


.B 


with  a  similar  expression  for  1  .  The  corresponding  limits  on  the 


-B 


y   integration   in  I     will  be   from  -y     to   0.    In  the  end  we  shall 


take  the  limits  y    =  0  and  y    =  oo. 

J  0  oo 
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Now  we  evaluate  the  q-integral  by  using  the  contour 
integration  as  in  the  case  of  integrals  in  eqs  (C.19)  and  (C.31). 
First,  consider  the  integral  of  the  lower  y  limit,  the  term 
containing  y  .  In  this  limit,  as  discussed  earlier  in  this 
appendix,  we  have  to  take  the  contour  in  the  LHP  as  shown  in 
figure  1-11.  Thus  we  find  that  the  integral  will  contain  only 
positive  powers  of  yQ  and  therefore  will  be  zero  in  the  limit  yQ  = 
0. 

Now  consider  the  integral  of  the  term  containing  y  .  In  this 
case  we  have  to  choose  the  contour  in  the  UHP.  We  have  to  include 
the  residue  at  Q  ,  the  zeroes  of  D(q),  and  at  qQ  =  l/2,  which 
arises  from  the  factor  (q-L/2)  in  the  denominator.  For  the  moment, 
ignore  the  pole  at  qQ  =  l/2  The  integral  will  then  contain  terms 
with  (t.Qr+0.5)  as  powers  of  y^.  In  this  case,  as  given  by  eqs 
(C.17)  and  (C.18),  Im(Qr)  ^  0.  Obviously  the  terms  for  which  Im 
(Qr)  -  0.5  will  also  vanish  in  the  limit  y^  =  ».  Hence  the  only 
terms  which  will  contribute  to  the  integral  are  those  for  which 


0  ^  Im(Qr)    <  0.5.  (C.35) 

Now  we  consider  the  special  term  for  which  q  =  l/2  for  which 
the  exponent  of  y^  is  0.  As  shown  in  section  3A,  the  value  of  D(q) 
is  also  zero  at  q  =  l/2.  However  it  does  not  represent  a 
second-order  pole  because,  as  given  by  eq  (97),  P(q)  is  also  zero 
at  q  =  l/2.  Therefore,  the  integrand  in  eq  (C.34)  has  only  a 
first-order  pole  at  qQ  =  l/2.  The  contribution  of  this  pole  to  the 
integral  can  be  written  as 

IAQ  =    2n  Lim       [  -4$-]-  <c-36> 
q=L/2    L  J 

Similarly,     we     obtain     the     following     expression     for  the 
integral  of  the  force  on  the  negative  Y  axis: 
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I 


B 


0 


=    2tt  Lim 

q=L/2 


(q) 


(C. 37) 


Adding  eqs  (C.36)  and  (C.37),  and  using  eqs  (93)  and  (94)  in  the 
definition  of  P(q)  as  given  by  eq  (55)  ,  we  find  that  the  numerator 
of  the  limit  will  approach  zero  faster  than  the  denominator.  The 
total  contribution  to  the  integral  from  the  pole  at  qQ  =  i/2  will, 
therefore,  be  zero. 

For    satisfying    eq    (12)  ,    the    only    contribution  to    1^.  can 

therefore  come  from  the  poles  at  q  =  Qr,   where  Qr  is  subject  to 

the  constraint  given  by  eq  (C.35).  Hence,  we  include  only  these 
contributions,  so 


with  a  similar  expression  for  I  .  The  total  integral  of  the  force 
function  is  given  by  eq  (C.33). 

If  a  Qr  does  not  exist  which  satisfies  (C.34),  then  1^  is 
zero  and  eq  (12)  is  satisfied.  If  1^  is  nonzero,  then  additional 
force  functions  have  to  be  applied  to  the  solid  corresponding  to 
the  same  Qr  which  give  a  finite  contribution  in  eq  (C.38)  so  that 
eq    (12)    is  satisfied.    These  terms  will  give   zero  at  x  =  0, 

since  they  are  the  solution  of  the  homogeneous  Hilbert  problem. 

The  expression  on  the  RHS  of  eq  (C.38)  diverges  for  y^  =  ». 
The  singularity  will  cancel  from  the  force  integral  calculated 
from  the  final  solution.  The  force  integral  will  be  zero,  since  we 
have  ensured  that  eq  (12)   is  satisfied  by  the  final  solution. 


Qr*L/2 


0.5 


(C. 38) 
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Free  Surface 


fA(y) 


/ 
/ 
/ 
/ 
/ 
/ 


Solid  A  (UHP) 


Interface 


/ 
/ 

fB(y)  / 

/ 
/ 


Solid  B  (LHP) 


Figure  1-1:     A  bimaterial  composite  containing  a  free  surface  normal  to  the 
interface  and  the  coordinate  system  used  in  these  calculations.     The  force 
functions  FA(y)  and  FB(y)  are  applied  just  outside  the  free  surface. 
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Figure  1-9:  The  contour 
function  HA(z) .  The  cut 
angle  is  0  on  the  branch 


used  for  evaluation  of  the 
PQRS  is  along  the  positive 
PQ  and  2tt  on  the  branch  RS 


integral  for  the 
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Figure  I- 10:     The  contour  used  for  evaluation  of  the  integral  for  the 
function  HB(z) .     The  cut  PQRS  is  along  the  negative  real  axis.     The  phasf 
angle  is  n  on  the  branch  PQ  and  -tt  on  the  branch  RS . 
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Figure  I -11:     The  semicircular  contours  used  for  evaluation  of  the 
q-integrals  in  the  UHP  (solid  line;   for  mod  z  >  1)  and  in  the  LHP  (dotted 
line;   for  mod  z  <  1). 
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Generalized  Plane-Strain  Analysis  of  a  Bimaterial  Composite 
Containing  a  Free  Surface  Normal  to  the  Interface 

by 

V.K.Tewary  and  R.D.  Kriz 
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1.  Introduction 


In  this  part,  we  extend  the  elastic  plane-strain  Green's 
function  calculated  in  Part  I  to  account  for  the  generalized  plane 
strain.  We  apply  it  to  calculate  the  stress  and  the  displacement 
field  in  a  bimaterial  composite  containing  a  free  surface  normal 
to  the  interface  which  is  subjected  to  generalized  plane  strain. 
The  generalized  plane  strain  is  induced  by  subjecting  the  solid  to 
an  out-of -plane  load  such  that  the  strain  component  c33  is 
constant. 

The  result  is  obtained  in  terms  of  a  closed  integral 
representation,  which  is  evaluated  analytically  as  well  as 
numerically.  As  in  Part  I,  the  analytical  evaluation  of  the 
integral  leads  to  a  series  representation — not  necessarily 
convergent — of  the  stress  from  which  the  singularities  in  the 
stress  can  be  identified  precisely.  The  numerical  calculation  of 
the  integral  representation  does  not  suffer  from  convergence 
problems.  It  contains  singular  as  well  as  nonsingular  terms  and 
gives  reliable  values  of  the  stress  and  the  displacement  field  in 
the  entire  solid.  The  method  is  applied  to  a  cubic  solid 
containing  a  Z-5  grain  boundary  and  to  fiber-reinforced  laminated 
composite  with  different  relative  fiber  orientations. 

In  general,  the  nature  of  singularities  in  the  generalized 
plane-strain  problem  is  similar  to  that  obtained  in  Part  I.  The 
power  of  a  singularity  depends  upon  the  elastic  constants  of  the 
constituent  materials,  whereas  its  weight  depends  upon  the  nature 
of  the  loading  as  well  as  the  elastic  constants.  In  some  cases  the 
singular  term  may  also  have  a  logarithmic  oscillatory  factor. 

The  plan  of  this  part  is  as  follows.  In  section  2,  we 
reformulate  the  Green's  function  method  for  generalized 
plane-strain  problems  and  apply  it  to  calculate  the  displacement 
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and  stresses  in  an  infinite  (without  a  free  surface)  bimaterial 
composite  subjected  to  an  out-of -plane  load.  The  displacement  and 
the  stress  components  t  .  thus  calculated  satisfy  the  required 
continuity  conditions  at  the  interface.  We  then  use  these  results 
and  the  technique  qiven  in  Part  I  to  calculate  the  stress  and  the 
displacement  field  in  a  similarly  loaded  composite  which  contains 
a  free  surface.  These  calculations  are  given  in  section  3.  The 
analysis  is  applied  to  a  Z-5  tilt  grain  boundary  in  a  cubic  solid 
in  section  4  and  to  laminated  composites  in  section  5.  A  summary 
and  a  discussion  of  the  results  are  presented  in  section  6. 

2.    Displacement   and   the    Stress   Field   in   an    Infinite  Composite 
Subjected  to  Generalized  Plane  Strain 

In  this  section  we  calculate  the  displacement  and  the  stress 
field  in  an  infinite  bimaterial  composite  which  has  no  free 
surfaces.  A  schematic  representation  of  the  model  solid  considered 
in  this  section  is  given  in  figure  II-l,  which  also  shows  the 
coordinate  axes.  This  model  is  same  as  used  in  Paper  I.  The 
notation  and  the  coordinate  axes  used  in  this  paper  are  the  same 
as  in  Part  I.  The  notation  in  this  paper  is  mostly  the  same  as 
used  in  Paper  I . 

The  displacements  in  both  UHP  and  LHP  have  to  satisfy  the 
following  equation  of  elastic  equilibrium  (add  superscripts  A  or  B 
corresponding  to  UHP  or  LHP  respectively) : 

a2  u. 

cikjl   3_    =  o.  (1) 

dXk  SX-j^ 

Equation     (1)      is     the     homogeneous     part     of     the  elastic 
equilibrium  equation  valid  in  the  absence  of  any  body  forces.  As 
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shown  in  Paper  I,  its  solution  in  a  particular  region  can  be 
written  using  the  elastic  Green's  function  as 

u(x)   =        G(x,x')f(x')   dx'  ,  (2) 
LJ 

where  f(x')  is  an  arbitrary  function  which  is  applied  in  a  region 
L  such  that  the  region  L  is  outside  the  region  of  solution  of  eq 
(1)  .  The  integration  in  eq  (2)  is  carried  out  over  the  entire 
region  L.  The  function  f (x' )  in  eq  (2)  is  usually  referred  to  as 
the  force  function.  The  force  function  in  eq  (2)  does  not 
represent  a  physical  force  on  the  solid.  It  is  a  mathematical 
artifact  which  is  used  to  satisfy  the  prescribed  boundary 
conditions . 

In  the  plane-strain  problem,  u(x),  G(x,x')  and  f(x')  are 
functions  of  only  the  x1~  and  x^-  coordinates  and  are  independent 
of  the  x.j-  or  the  Z-  coordinate.  Even  in  plane-strain  problems, 
the  Z-  components  of  these  quantities,  in  general,  will  not  be 
zero.  The  solution  of  eq  (1)  using  eq  (2)  for  the  plane-strain 
case  has  been  discussed  in  Paper  I. 

In  the  generalized  plane-strain  problem,  the  solid  is 
subjected  to  a  load  in  the  Z-  direction  such  that  the  e33  is 
constant  everywhere  in  the  solid  and  has  a  prescribed  value,  e 
being  the  strain  tensor.  The  displacement  field  in  the  solid  can 
therefore  be  written  in  the  form 

Vl^)   =  ^  *l3  x3  +  up.(x),  (3) 
where  t)  is  the  prescribed  value  of  e33  so 

H  =  c33  (4) 
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is  constant  and  u  (x)    is  a  vector  which  is  a  function  of  only  x„ 

-p  -  J  1 

and  x2  and  not  of  x3 .  The  subscript  «  on  u  indicates  that  it 
denotes  the  displacement  field  in  an  infinite  composite  with  no 
free  surfaces.  Similarly,  the  subscript  p  on  u  on  the  RHS 
indicates  that  it  is  a  solution  of  the  plane-strain  problem  as 
described  below. 

Following  the  notation  used  in  Part  I,  x,   x'  will  denote  a  2D 

(2  dimensional)   vector  on  the  XY-  plane  with  Cartesian  components 

x^f    x2   and  x'  ^,    x' 2  respectively.   The  3D  vectors  will  be  denoted 

by    the    upper    case    variables    X,    X'  .    For    example,    the  position 

vector  X  in  eq   (3)   denotes  a  3D  vector  with  Cartesian  components 

x,  ,    x0   and  x_ .   Our  object  in  this  section  is  to  determine  u  (x) 
1'      2  3  J  -pv-' 

such  that  um(x)  is  a  solution  of  eq  (1)  and  it  satisfies  the 
following  boundary  conditions, 

u*.    [x1   ,0]   =  u*.    [x1   ,0]        (-co  *  Xl  s  .)  (5) 

and 

T«2i  rxi  '°]  =  T!2i  txi  '°3     ("M  *  xi  *  <6> 

where  Tm2  is  a  vector  which  is  defined  in  terms  of  the  components 
t .     of  the  stress  tensor  (for  superscripts  A  and  B)  as 


T»2i  =  Ti2    •  <7' 


The  stress  components  can  be  calculated  from  the  displacement 
field  by  using  the  formulas, 


Su,  flu.  5u, 

Til  =  Cillk  —  +  Cil2k  —  +  Cil3k  —     '  <8> 
SX1  9X2  aX3 
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Suk  Suk  8uk 

Ti2  =  Ci21k  ^  +  °i22k     +  Ci23k  173  ■  (9) 

and 

Ti3  =  Ci31k  —  +  Ci32k  —  +  Ci33k  —  ' 

9X1  ax2  5X3 

The  first  term  in  eq  (3)  is  a  constant  with  respect  to  the 
differential  operators  in  eq  (1) .  It  would  therefore  automatically 
satisfy  eq  (1).  Since  u  (x)  depends  only  on  x  and  x  ,  we  can  use 
the  plane-strain  Green's  function  obtained  in  Paper  I  to  determine 
this  part  of  the  displacement  field. 

First  we  consider  the  displacement   field   in  the  solid  A  in 

the    UHP.     Following    the    method    given    in    Paper    I,     we  apply 

hypothetical  forces  just  outside  the  region  A  at  the  set  of  points 

(x1#   -e)  where  x^j^  varies  from  -co  to  +  oo  and  c  is  a  small  positive 

constant  which  will  be  taken  to  be  zero  in  the  limit.  The  positive 

value  of  e  ensures  that  the  hypothetical  forces  are  applied  just 

outside  the  region  of  solution.   We  denote  this  force  function  by 
A 

f  (x1)  .  By  using  eq  (2)  and  as  given  in  Paper  I,  we  can  write 
Up(x)  in  the  form, 


,00 


HpA(£)  =  GA(x;x^-c)fA(x^)   dx^  ,  (11) 

—CO 

where  G  (x;x')  denotes  the  Green's  function  for  the  solid  A  when 
the  solid  B  is  not  attached  to  it. 


As  described  in  the  previous  paragraph,  the  Y  component  of  x' 
has  been  taken  to  be  -e  to  ensure  that  the  force  f  (x^)  has  been 
applied  just  outside  the  region  A.  Equation  (11)  then  gives  a 
solution  of  eq  (1)   in  the  entire  space  of  A. 
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The  Green's  function  for  the  solid  A,   in  the  absence  of  solid 
B,  can  be  written  as  follows(  see,  for  example,  Paper  I): 

GA(x;x')  «  -  JL  £  yA(a)   ln(zA  -  pA)  ,  (12) 


where 

za  -  xl  +  P«  x2  '  <13> 

P*  "  Xi  +  P«  X2    •  <14> 

A  A 

and    r   (a)    and    pa    have  been    defined    by    eqs    (A.  17)    and    (A.  19) 

respectively  in  Appendix  I-A  of  Part  I.  From  eqs  (3)  and  (12)  we 
obtain 


uA(X)   =  uQ     -  -A-  £  yA(a)  ln(zA  -  t-Le) f A(t) dt,  (15) 

a  -coJ 


where  it  is  understood  that  only  the  real  part  of  the  RHS  in 
eq  (15)  represents  the  displacement  field  and  where 

U0i  =  77  5i3   X3    *  (16) 


In  the  arqument  of  the  loqarithmic  function  on  the  RHS  of  eq 

A 

(15)  we  have  written  p    e  =  te.  This  is  because  we  have  to  put  e  = 

•  A  • 

0     in    the    limit    and    therefore    the    maqnitude    of    p^     is  not 

siqnificant.  The  only  thinq  which  could  possibly  be  siqnificant  in 

the  loqarithmic  term,    as  in  Part  I,    is  the  siqn  of  the  imaqinary 

part  of  pa  which,   as  qiven  in  Appendix  I-A,    is  positive.  However, 

in  the  present   case,    we   do   not   need   to    introduce   a   cut    in  the 

complex  plane  so  that  the  limit  e  =  0  can  be  taken  immediately. 
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Following  the  method  given  in  Papers  I  and  II  and  also  Part 
I,    and  using  eq    (9),   we  obtain  the  following  expression  for  the 


stress 


£2<*>   -  202  -  '4-  I  2   <«>      J         A       ~       '  <17> 

a  -ooJ  Z       -  t+LC 

a 


where  the  vector  TQ  is  defined  as 


T02i  "  Ci233  *   •  <18> 


In  a  similar  fashion  we  obtain  the  following  expressions  for 
isplac 
in  the  LHP: 


the  displacement  field  and  the  stress  component  t . 9  in  the  solid  B 


00 

(X)  =  uQ     -  ][  iB(oc)  ln(z*  -  t+te)fB(t)dt,  (19) 

a  -ooJ 


fB(t)dt 


=^02  -^rl?w   I    -b    -    >  <20> 

a  -ooJ         z      -  t+LC 

a 


where 


S02i  "  C!233   *   •  <21> 
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Now  we  apply  the  boundary  condition  given  by  eq  (5)  on  eqs 
(15)  and  (19)  .  This  gives  the  following  relation  between  the  two 
force  functions: 


Is  £A(t)  =  *s  £B(t)* 


(22) 


A  B 

For  convenience,  we  write  f  (t)  and  f  (t)  in  terms  of  another 
force  function  fQ(t)  as 


^(t)  =         [Ig]"1  fQ(t) 


(23) 


and 


fB(t)  =  cig]"1  f0(t) 


(24) 


A  B 

We  now  substitute  for  f  (t)   and  f  (t)    from  eqs   (23)   and  (24) 
in  eqs   (17)   and   (21)   and  apply  the  boundary  condition  given  by  eq 
(6) .  After  rearranging  the  terms,  we  obtain  the  integral  equation 


TT 


-co 


xi  "  t 


=  5  t0  , 


(25) 


where 


D  = 


-B     -B  -1 


A   .   A  -1 

a     [If  ] 


-1 


(26) 


T     =  TB 
-0  -02 


TA 

-02  ' 


(27) 
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and  the  integral  on  the  LHS  of  eq  (25)  is  a  principal  value 
integral . 


Equation  (25)   is  in  the  standard  form  of  a  Hilbert  transform. 
Its  solution  is  given  by 

-00 


dt 

  ,  (28) 


—00 


xl  - 1 


where  we  have  used  the  fact  that  the  RHS  of  eq  (25)  is  a  constant. 
The  integral  on  the  RHS  of  eq  (28)  can  be  easily  carried  out.  For 
this  purpose,  it  would  be  convenient  to  evaluate  the  integral 
taking  the  upper  and  the  lower  limits  to  be  t  and  -t 
respectively  and  then  taking  the  limit  t    -  oo.  The  result  is 


fQ(x1)   =     i  D  Tq   ,  (29) 

where  we  have  taken  the  following  limiting  values  of  the 
logarithmic  function: 


ln(x.   +  t  )   =  ln(t    )  (30) 

x      1  00  x     03     '  x  ' 


and 


ln(x,   -  t  )  =  ln(  t    )   +  til  (31) 

x    1  oo'  x       oo    '  x/ 


in  the  limit  t    -»  oo. 

00 

From    eqs     (23)  ,     (24)     and     (29)  ,     we    obtain    the  following 
solution  for  the  two  force  functions: 

fA(t)   =  -L  M  TQ  (32) 
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and 


fB(t)  = 


(33) 


where  the  matrices  M  and  N  have  been  defined  in  eqs  (A.  11)  and 
(A. 12)   in  Part  I. 

We  see  from  eqs  (32)  and  (33)  that  the  two  force  functions 
are  constants.  Using  these  results  in  eqs  (15)  and  (19)  gives  the 
required  expressions  for  the  displacement  field  in  the  solid.  The 
logarithmic  integrals  in  eqs  (15)  and  (19)  can  be  easily 
evaluated.  Thus  we  obtain  the  following  result  for  the 
displacement  field  in  the  UHP  and  the  LHP: 

uA(X)   =  u0  -  Y  «rA(a)M  Tn  zA 

-oo v -0      L  -  v   '  -  -0  a 

a 

and 

uB(X)   =  un  -  V  yB(a)N  T„  zB  , 

-oo x           -0      L  -   v   '  -  -0  a 

a 


(34) 


(35) 


where  we  have  used  eqs  (30)  and  (31)  and,  as  is  usual  in  continuum 
mechanics,  ignored  the  terms  which  correspond  to  the  rigid  body 
displacements  in  the  solid. 

For  later  use,   the  expressions  for  the  stress  components 
and  whicn  can  be  obtained  from  eqs   (34)   and  (35)   by  using  eqs 

(8)   and  (9) ,   are  given  here. 

ti  -  loi  -      3  To  '  (36> 


TB  =  TB  -  QB  N  T  (37} 
-ool       -01       -s  -     0    '  {  ' 
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and 


^2  "  ^02   "  4  S  T0  <38> 


£2  m  So2  "  5  H  To  '  <39> 


A  B  A  B 

where  the  vectors  and  are  defined   in  analogy  with  eqs 

(7)   and  (18),  as 

£li«>  =  TiiB  <40> 

and 

mA ,  B        A ,  B  ,  .  _  . 

2dii  =  cil33  71     •  <41> 

Equations   (34),    (35),    (38)   and   (39)   satisfy  the  two  boundary 

conditions  given  by  eqs   (5)   and   (6)  .  Moreover,   the  second  term  in 

both  the  eqs    (34)    and   (35)    is  a  function  of  x1  and         and  not  of 

x3 ,    whereas   the   first  term   is   a   function   of   x3   only.    Hence  the 

value  of  the  strain  component  e33  will  be  17  as  prescribed  by  eq 

(4) .   Equations   (34)   and   (35)   are,   therefore,   the  desired  solution 

for  the  displacement   in  an   infinite   composite   solid.    It  may  be 

noted  that  x . _  is  constant  in  the  solid. 
12 

3.  Displacement  and  the  Stress  Field  in  a  Composite  Solid 
Containing  a  Free  Surface 


In  this  section  we  shall  apply  the  results  obtained  in  the 
last  section  to  calculate  the  displacement  field  and  the  stress 
distribution  in  composite  containing  a  free  surface.  The  model 
solid  and  the  coordinate  axes  are  shown  in  figure  1-1  of  Part  I. 
Mathematically,  the  difference  between  the  calculations  of  Part  I 
and  the  present  calculations  arises  because  in  Part  I  a  line  load 
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was  specified  at  a  given  point  whereas  in  the  present  case  the 
strain  e33  is  specified. 

We  start  by  writing  the  displacement  fields  in  A  and  B  in  the 
following  form,  which  is  analogous  to  that  given  in  eq  (3) . 


uA(X)   =  uA.(X)   +  uA.(x)  (42) 

lv-'  oolv-7  Slv-'  v  ' 

and 

uB(X)  =  uB.(X)   +  uB.(x)    ,  (43) 

1    —  col  v  -'  SI v -'     '  v/ 


A  B 

where  u  '    (X)    are  the  displacement   field  vectors    in   an  infinite 

composite  with  no  free  surfaces  as  given  by  eqs   (25)   and  (36)  ,  and 
A  B 

u  '    (x)   are  functions  of  only  x.   and  x_  ,  and  not  of  x_ . 
-s      —  12  3 

In    view    of    the    remark    made    at    the    end    of    the  preceding 

section,    the  strain  tensor  calculated  from  displacement  field  as 

given    by    eqs     (42)     and     (43)     will    satisfy    eq     (4) .     Our  task, 

A  B 

therefore,   is  to  determine  u  '    (x)   such  that  eqs   (42)   and  (43)  are 

™*  s  ^™ 

solutions  of  the  elastic  equilibrium  equations  [eq  (1)  with 
appropriate  superscripts]  in  the  HHP  and  LHP  respectively  and 
satisfy  the  boundary  conditions, 


uAj,[x1,0]  =  uB.[x1,0]  (-co  <  Xl  *  »),  (44) 
xAi2[x1,0]   =  zBLl   [x1,0]        (-co  <  x±  *  co),  (45) 


and 


T  ilC°'X2]   =  °  (°  "  X2  "  M)  '  (46) 


TB±1[0,x2]   =0  (0  *  x2  ^  -oo).  (47) 
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We  can  determine  u  '  (x)  by  using  the  plane-strain  Green's 
function  as  quoted  in  Appendix  I-A.  This  Green's  function  is  a 
solution  of  eq  (1)  and  satisfies  the  boundary  conditions  given  by 
eqs  (44)  and  (45) .  In  order  to  satisfy  the  boundary  condition 
given  by  eqs  (46)  and  (47)  ,  we  shall  use  the  same  technique  as 
used  in  Part  I . 

3A.   Solution  in  Terms  of  an  Integral  Representation 

As   in  Part  I,   we  apply  a  hypothetical  distribution  of  line 

A 

forces  just  outside  the  free  surface  denoted  by  F  (t)  in  the  UHP 
and  F  (t)  in  the  LHP.  As  shown  in  figure  II-2,  these  forces  are 
applied  just  outside  the  free  surface  at  continuous  set  of  points 
x1=-e  and  x2  =t,  where  e  is  a  small  positive  number  which  is  taken 
to  be  zero  in  the  limit.  The  range  of  the  real  variable  t  is 
(0  *  t  *  oo)   for  FA(t)   and  (-«  <  t  ^  0)   for  FB(t)    in  the  LHP. 

Thus,    using   the   definition   of   the   Green's    function   and  as 

given    by    eq    (2)  ,    we    can    write    the    following    expressions  for 

uA'B(x)   in  UHP  and  LHP 
s  - 

oo  0 

Ug.(x)  =     I  GQAA(x;t)FA(t)dt  +        GQAB(x;t)FB(t)dt  (48) 
0 

and 

CO  0 

uB.(x)   =      '  GQBA(x;t)FA(t)dt  +        GQBB(x;t)FB(t)dt/  (49) 

~  00 

0 

where  GQ  denotes  the  Green's  function  for  an  infinite  bimaterial 
composite  which  is  given  by  eqs  (A.l)  —  (A.  4)  in  Appendix  I-A. 
Using  these  equations  and  eqs  (117)  and  (118)  of  Part  I,  we  obtain 
the  following  expressions  for  the  displacement  field  from  eqs 
(42) ,    (43) ,    (48)   and  (49) : 
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ua(x)  =  u^(x)  +  l 


IA(a) 


r    Z    +   C  ^ 

a 


zA+  e 


a 


(X 


+  L 


-A 

z 


Sp  L  1  Pp  J  1  Pp  J 


a0  l 


*AB<a'^  5 


B 


a 


K  P 


B 
$ 


+  EAB(a#/3)  U 


B 


(50) 


and 


(for  0  ^  x2  ^  oo) 


u  (x) 


-B 

z 


=  H*(X)   +  t  £    EBA(a^)  UA[  ^]   +         ( a ,  /3 ) 
a/3  L 


B 

f  z~  \ 


a 


K  P 


+  L 


a  u 


B  /    x  TTB 

y  (a)  y 


'  za+  g  - 

B 


+  rB(a)  uB 


-B , 
Za+  C 
-B 


-  L 


a/3  L 


EBB(a,P)  U 


B 


-B 

Z  _.  \ 


a 


v  P 


B 
0 


+  EBB(a,0)  H 


B 


B 

za 


(51) 


(for  -oo  <  x2  ^  0)  , 

where,  as  in  Part  I,  the  U-  functions  are  defined  as  the 
indefinite  integrals, 


U^z)  =      HA(z)dz  , 


(52) 
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UB(z)   =      HB(z)dz  , 


(53) 


HA(z)  = 


FA(t) 


27TL 


0J 


and 


HB(Z)  = 


t  -  z 


-°  FB(t) 


dt 


27TL 


dt, 


(54) 


(55) 


t  -  z 

and  the  other  symbols  have  been  defined  in  Part  I. 

For  later  use,  we  also  write  below  the  stress  components 

' ' 1'  w^^c^  can  ke  easi] 
given  in  eqs   (50)   and  (51) . 


and  t • «i  which  can  be  easily  obtained  from  the  displacement  field 


To(x)  =  TA 


2oo 


+  I 


]T    {**(«)/  pa}  5A 


a 


f  Z"  1 

A 

v       P  * 

+  j2A(a)/  pA|  HA 


-A . 
Zcc+  c 
-A 


a/3 


DAA(a,<3)  H 


A 


-A 


-A 

z 

a 

^  PA  ' 


-  L 


a/3  l 


DAB(a,P)  H 


B 


A 


'  P, 


B 


+  DAB(a,f3)  H 


B 


-A 

z 

a 


(56) 
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(x)   =  TB 


2  00 


+  L 


+  L 


r  -B  B 

L    eBA<«^)  sA(  ^-)  +  eBA(a,P)  „M  !«_] 


I  [       p£}  sB(       ♦  {?(->/  PaB}  bb(  !llLj 

Pa/ 


a 


-  L 


r  -B  B 


(57) 


ool 


«  L  Pa    ^  1      Pa  J 


+  i 


a/3 


.A'  2 


A 
a 


Saa'"'*3)   !T     4~     +  K      (a,P)   HAf  -I-' 
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-  L 


a|3  L 


KBB(«,0)  H 


B 


-B 


k  P 


B 

0 


+  KBB(a,/3)  H 


B 

B(  fa_  1 
-B 

1  P0 


(59) 


where  the  vectors  and  Tm2   are  given  by  eqs    (3  6) -(39)    and,  in 

an  analogous  manner,   we  have  defined  the  vectors  T. (x)    and  T_  (x) 

th  — 
such  that  their  i      components  are  equal  to,   respectively,  z±i^^ 

and  ri2 (x) . 

As  mentioned  earlier,    eqs    (50)    and    (51)    will  satisfy  eq  (1) 

everywhere    in    the    region    of    solution    and    also    the  continuity 

conditions  at  the  interface  as  given  by  eqs   (44)    and   (45)    for  all 
A  B 

values  of  F  (t)  and  F  (t)  .  This  is  because  the  forces  have  been 
applied  outside  the  region  of  solution  and  the  Green's  function 
used  in  writing  eqs  (50)  and  (51)  satisfies  eqs  (44)  and  (45) .  As 
in  Part  I,  we  have  to  determine  these  two  force  functions  so  that 
the  free  surface  boundary  conditions  as  given  by  eqs  (46)  and  (47) 
are  also  satisfied.  In  addition,  the  force  functions  have  to 
satisfy  the  following  condition  required  for  the  displacements  to 
be  single-valued: 


00  (J 

FA(t)   dt  +  FB(t)   dt  =  0 


(60) 


—  00 


The  H-  functions  are  determined  in  the  same  way  as  in  Part  I. 
We  write  these  in  the  form, 

HA(z)  = 

—00 


TTA,   .    iq-0.5  , 
\f  (q)  z  dq 
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and 


YB(q)z^-°-5  dq     ,  (62) 

—00 

where  (q)  and  VB(q)  are  the  unknown  functions  of  the  real 
variable  q,  which  have  to  be  determined. 

Following  exactly  the  same  procedure  as  given  in  Part  I  and 
using  the  same  block  representation,  .  we  obtain  the  following 
result  for  the  vectors  V^q)   and  VB(q), 


HB(z)  = 


M(q)   V(q)   =  N(q) 


(63) 


where  M(q)  is  exactly  the  same  as  given  by  eqs  (38) -(41)  of  Part  I 
and  the  block  elements  of  N(q)  are 


Na(q)  =  - 


and 


NB(q>  =  - 


-03  1 

2m 


-col 
2TTL 


00 


-Lq-0.5  , 

y   H  dy 


(64) 


0J 


-00 


I y I     *         dy . 


(65) 


The  integrals  in  eqs  (64)  and  (65)  can  be  evaluated  in  an 
elementary  manner.  However,  in  order  to  be  able  to  handle  the 
apparent  divergence  at  oo,  we  evaluate  the  integral  in  eq  (64) 
between  the  limits  0  and  y^  and  that  in  eq  (65)  between  -y^  and  0, 
and  in  the  end  take  the  limit  y  =  «.  Thus  we  obtain  the  following 
expressions  for  the  block  elements  of  N: 
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NA(q)   =  t£(q)7(-Lq+0.5) 


(66) 


and 


NB(q)   =  t®(q)/(-Lq+0.5) , 


(67) 


where 


and 


tA(q)  = 

—  oo  -1 


tB(q)  = 

—  00  x 


-oo  1  -Lq+0.5 

^00 

2TTL 


-oo  1  -Lq+0.5 

^00 

27TL 


(68) 


(69) 


The  transform  of  a  constant  in  q-space  as  given  by  eqs  (68) 
and  (69)  is  singular  in  the  limit  =  oo.  However,  as  shown  in 
Appendix  II-A,  this  singularity  cancels  out  when  the  inverse 
transform  is  taken. 

The  solution  of  the  matrix  equation  (63)  can  be  easily 
written  as 


Y(q)   =  [M(q)]   1  tM(q)/(-Lq+0.5) . 


(70) 


For  later  use,  we  write  the  following  block  structure  of  eqs  (63) 
and  (70)   as  in  Part  I: 


-AA 


M 

^  -BA 


M 

-AB 


M 

-BB 


1 

Y*(q)' 

1 

J 

,  YB(q). 

(-Lq+0. 5) 

A  1 


—  oo 
.B 

•oo 


(71) 
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and 


I  YB(q)J 


(-Lq+0.5)  D(q) 


C 

v  -BA 


r  a  i 

tA 

Sab 

—  00 

c 

tB 

-BB  J 

<     -00  ' 

(72) 


where  C  is  the  matrix  of  cofactors  of  M. 

As  remarked  in  Part  I,  the  representation  of  V(q)  in  terms  of 
the  cof actor  matrix  as  given  by  eq  (72)  is  useful  mainly  for  the 
purpose  of  a  qualitative  analysis  of  the  nature  of  H(z)  .  For 
numerical  calculations,  it  is  better  to  use  eq  (70)  ,  which 
involves  only  the  inversion  of  a  6  x  6  matrix. 

Equation  (70)  gives  the  particular  solution  of  the  Hilbert's 
equation.  The  general  solution  is  obtained  by  adding  the 
homogeneous    and    the    particular    solutions    and    is    written  as 


H(z)  = 


(     A  } 

3p<z> 

+ 

r  a  i 

.  3o<z>. 

(73) 


where  the  subscript  p  denotes  the  particular  solution  which  is 
given  by  eqs  (61)  and  (62)  along  with  eq  (70)  .  The  homogeneous 
part  of  the  solution  in  eq  (73) ,   labelled  with  the  subscript  0,  is 


H0(z)   =£h(Qr)   zLQr-°-5  , 


(74) 


where  the  vector  h  is  a  solution  of  the  matrix  equation, 


M(Qr)h(Qr)  =  0 


(75) 
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and  Qr  is  a  root  of  the  determinantal  equation, 

D(Qr)   =  0   .  (76) 

One  of  the  components  of  the  vector  h  for  each  Q  will  be 
arbitrary.  It  has  to  be  chosen  so  that  eq  (60)  is  satisfied.  These 
constants  can  also  be  chosen  to  satisfy  any  other  boundary 
condition  on  the  stress  field,  such  as  those  which  may  be  required 
at  the  outer  surfaces  in  a  finite  solid. 

3B.  Series  Representation  for  the  Particular  Solution 

In  this  subsection,  we  shall  obtain  the  series  representation 
for  the  particular  solution  for  H(z) .  It  may  be  remarked  that  for 
an  actual  calculation  of  H(z)  ,  it  would  be  generally  more 
convenient  to  evaluate  numerically  the  integral  over  q  in  eqs  (61) 
and  (62).  The  series  representation  for  H(z)  may  not  give 
numerically  accurate  values  of  H(z)  except  for  limiting  values  of 
mod(z),  because  the  convergence  of  the  series  is  not  guaranteed. 

However,  the  series  representation  for  H(z)  serves  two  useful 
purposes:  (i)  it  can  be  used  as  a  numerical  check  on  the  direct 
calculation  of  H(z)  obtained  by  numerical  integration  in  eqs  (61) 
and  (62)  ,  and  (ii)  it  helps  in  identifying  the  singularities  in 
H(z)  at  mod(z)  =  0,  which  provides  a  better  physical  understanding 
of  the  nature  of  the  stress  distribution  in  the  solid.  Moreover, 
as  mentioned  in  section  1,  a  knowledge  of  the  singularities  in 
H(z)  is  essential  for  a  proper  numerical  modelling  of  the  solid  by 
using  finite  element  or  other  such  techniques. 

First,  we  shall  consider  H  (z)  .  We  see  from  eqs  (61)  and 
(72)   that  the  particular  solution  for  H  (z)   can  be  written  as 
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?A(q)  [z/yj 


iq-0.5 


dq 


(77) 


00 


(-Lq+0.5)  D(q) 


where 


PA(q)  = 


AB 


(78) 


with  a  similar  expression  for  H  (z) . 


We  evaluate  the  integral  in  eq  (77)  by  following  the  same 
procedure  as  used  in  section  3  in  Part  I.  Comparing  eq  (77)  with 
eq  (C.16)  of  Part  I,  we  find  that  the  only  poles  which  can 
contribute  to  the  integral  in  eq  (77)  are  those  arising  from  the 
zeroes  of  D(q)  and  the  single  pole  at  q__  =  -i/2.  The  set-I  poles 
in  eq  (C.16)  of  Part  I  arising  from  the  function  E(q)  are  not 
present  in  eq  (77)  except  for  the  single  pole  at  q  =  CL-^* 
Moreover,  in  the  two  examples  discussed  in  section  4  and  section 
5,  we  shall  see  that  D(q_1)  is  also  zero.  Hence,  the  pole  at  q  = 
is  a  second-order  pole  in  the  lower  half  of  the  complex  plane. 

Since  y^  approaches  oo  in  the  limit,  we  note  that  mod(z/yw)  in 
eq  (77)  will  approach  0.  Thus,  as  indicated  in  section  3B  and 
Appendix  I-C  in  Part  I,  we  have  to  choose  the  contour  in  the  LHP 
for  evaluating  the  integral  over  q  in  eq  (77) .  Thus,  by  using  eqs 
(C.34)  and  (73)  of  Part  I,  we  obtain  the  following  series 
representation  for  H(z) : 


lQ  -0.5+c 
r 


-  RH(z) , 


(79) 


r 


where  Q    is  a  root  of  eq  (7  6)  with  the  constraint, 
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Im  (Qr)   <  0  , 


(80) 


and  the  other  quantities  in  eq  (79)  have  been  defined  in  Part  I. 
In  eq  (79)  and  in  what  follows  in  this  section,  z  is  expressed  in 
units    of    y      for    notational    convenience.     In    actual  numerical 

00 

calculations  y    serves  as  a  scale  factor. 


00 


In  analogy  with  eq   (103)    of  Part  I,   we  write  eq   (79)    in  the 

form, 

(k  -0.5) 

Hp(z)  =  1  ^p(Qr)    [z]     r  exp[Lgrln(z) ]  +  RR(z) ,  (81) 


where  gr  and  -kr  are  the  real  and  imaginary  parts  of  Q^., 

Qr  =  gr  -  *.kr  (82) 

and 


Yp(Qr)  =  -  *(Qr)  g(i<Qr)   •  (83) 


We  see  from  eq  (81)  that  H(z)  will  have  the  same  kind  of 
singularity  and  the  logarithmic  oscillations  as  discussed  in 
section  3B  of  Part  I.  The  singularity  will  exist  only  if  a  Qr 
exists  such  that  its  imaginary  part  satisfies  the  constraint  that 


0  <  k     <  0.5   .  (84) 
r 

Similarly,  the  logarithmic  oscillatory  factor  will  exist  only  if 
the  real  part  of  Qr  is  nonzero. 

In    eq    (81)  ,    the    first    term    contains    sum    over    only  those 
values  of  Q    which  are  not  higher  order  poles.  The  contribution  of 
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higher  order  poles  is  included  in  Rjj(z) .  In  particular,  the 
contribution  R^Q(z)   of  the  second-order  pole  at 


q_i  =  Q_i  =  -l/2 


(85) 


to  Rjj(z)  ^as  been  given  below.  Following  the  method  leading  to  eq 
(104)    in  Part  I,  we  obtain 


d'  (q_x) 


E' (q_x)  +  ln(z) 


(86) 


Equation  (86)  shows  the  ln(z)  behavior  of  H(z)  in  agreement 
with  some  of  the  earlier  work  (see,  for  example,  [10]  and  other 
references  given  there) .  Further,  if  the  root  Qr  is  degenerate,  we 
shall  get  logarithmic  terms  with  higher  powers.  For  example,  if 
for  some  Qr  *  -l/2  which  satisfies  the  constraint  given  by  eq 
(84),  D(Qr)  has  a  factor  (q  -  Q^.)^ ,  then  by  using  the  Cauchy's 
theorem  we  can  easily  show  that  R„(z)  will  have  the  form, 


Rjj(z)  «  [z] 


(k  -0.5) 


„  ~  1 

[ln(z)]^"  exp[Lgr7?ln(z)  ] , 


(87) 


where        is  written  in  terms  of  its  real  and  imaginary  parts  as 


Q    =  g      -  t.k 


(88) 


with  the  constraint 


k      >  0  . 


(89) 


The  subscript  tj  in  eq  (88)  identifies  the  particular  value  of 
Qr  as  a  7] — fold  degenerate  root. 
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3C.  Force  Functions  and  the  Homogeneous  Solution 


In  this  subsection  we  calculate  the  force  functions  and  their 
integral  by  using  the  same  technique  as  in  Part  I.  By  using  the 
Plemelj  relation  as  given  by  eq  (B.42)  of  Part  I,  we  obtain  the 
following  expression  for  the  two  force  functions: 


ip(y)  = 


E(q)PA(q)dq 
(iq-0.5)  D(q) 


iq-0 . 5 


—  00 


(90) 


and 


i?(y) 


exp(irq)E(q)PB(q)dq 
(iq-0.5)  D(q) 


'y/y. 


iq-0 . 5 


(91) 


—  00 


where  P  (q)  has  been  defined  in  eq  (78)  ,  and 


PB(q)  = 


(92) 


The  integrals  over  q  in  eqs  (89)  and  (90)  can  be  evaluated  by 
using  the  same  method  as  given  in  Appendix  I-C  of  Part  I  for 
calculating  the  integral  in  eq  (C.31).  Unlike  eq  (C.31)  in  Part  I, 
the  integrands  in  eqs  (89)  and  (90)  do  not  have  the  factor  E(q)  in 
the  denominator  but  have  the  factor  (-iq  +  0.5).  The  singularity 
arising  from  this  factor  at  q  =  -t/2  cancels  out  with  the  zero  in 
E(q)  in  the  numerator.  The  integrands,  therefore,  have  only  a 
simple  pole  at  q  =  -l/2  which  is  a  zero  of  D(q) . 

We  need  only  to  calculate  the  total  integral  of  the  force 
functions  in  order  to  satisfy  eq   (60) .   The  total  integral  of  the 
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force  functions  has  also  been  evaluated  in  Appendix  I-C  of  Part  I. 
From  eq  (C.33)  of  Part  I  and  eqs  (88)  and  (89),  we  obtain  the 
following  expression  for  the  integral  of  the  force  in  UHP: 


Y 

CO 


ip(y)  dy 


,0°  Yn  E(q)  PA(q)dq 
.      D(q) (iq+0.5) (-iq+0.5) 

—  CO 

We    take    the    limit    as    Y     >    y     ->    co.    As    explained    in  the 

CO  CO 

Appendix  I-C  of  Part  I,  only  the  poles  in  the  LHP  contribute  to 
the  q-integral  at  the  lower  limit  of  y,  and  their  contribution  is 
zero.  At  the  upper  limit  of  y,  the  q-integral  is  shown  in  eq  (93) . 
In  this   case  the  poles   in  the  UHP  will   contribute   in  the  limit 

Y  =00. 
CO 

First  we  consider  the  pole  at  q  =  l/2.  In  contrast  to  eq 
(C.33)  of  Part  I,  the  integrand  in  eq  (93)  has  the  E(q)  term  in 
the  numerator,  which  cancels  the  singularity  arising  from  the 
factor  (-Lq+0.5)  in  the  denominator.  However,  as  given  in  eq 
(C.25)  of  Part  I,  D(q)  is  also  zero  at  q  =  l/2.  Thus  the  integrand 
in  eq  (93)  has  a  simple  pole  at  q  =  l/2  which  arises  from  the  zero 
of  D(q).  The  contribution  of  this  pole  is  given  by 

JT0     *  *n\  /  D'<V    '  <94> 

where  the  prime  denotes  differentiation  with  respect  to  q  and 

q0  =  l/2   .  (95) 


(93) 
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Similarly  we  obtain 


Using  eqs  (78),  (92),  (94),  (96)  and  eqs  (93)  and  (94)  of  Part  I, 
we  see  that  the  contribution  of  the  pole  at  qQ  to  the  total  force 
integral  is  zero. 

Further,  as  in  Part  I,  we  note  from  eq  (93)  that  the 
contribution  of  a  pole  for  which  ImfQ^J  >  0.5,  will  be  zero  in  the 
limit  =  oo,  because  it  would  lead  to  a  negative  power  of  Y^. 
Thus,  we  see  that  only  those  poles  will  contribute  to  the  integral 
which  satisfy  the  following  constraint: 


0  <  Im(Qr)   <  0.5   .  (97) 
This  leads  to  the  following  for  the  total  force  integral: 


lQ  +  0.5 

i  =  _  x      E(Qr)yw(YM  A.  )  r 


:,    =   2TTL    )  — 
J-^  D' 


D'  <V   <Qr  +  x/4) 

r  1 

[  PA(Qr)   -  Lexp(TTQr)   PB(Qr)   1,  (98) 

where  Q  is  subject  to  the  constraint  given  by  eq  (97)  .  If  a  Qr 
which  satisfies  this  constraint  does  not  exist,  then  1^.  is  0. 

As  remarked  in  Part  I,  the  apparent  singularities  in  1^  are 
not  relevant.  When  1^.  is  not  zero,  we  have  to  include  additional 
force  terms  corresponding  to  the  solution  of  the  homogeneous 
equation.  The  H-f unctions  as  obtained  by  these  terms  are 
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=    £  v£(Qr)zLQr-°-5 


(99) 


and 


H*(z)  =    £        V*  (Qr)   zLQr"  °'5  f  (100) 


where 


Qr*t-/2 


3£(QJ    =  2TTL 


0(Qr>    =  2TTL    (-LQ  4-0.5)0-  (Q   )  <101> 


and 


b  EB(Qr) 

W    -  2l"    (-.Qr+0.5)D'(Qr)    '  <102> 

If  we  separate  Q     into  its  real  and  imaginary  parts  as 

Qr0  =  ^ro  +  Lkr0  '  <103> 

then  we  see  that  the  H — functions  will  have  the  following 
dependence  on  z : 

-(k  +0.5) 

H(z)   *  (z)       ru          exp[cgr0ln(z)].  (104) 


3D.  Singularities  in  H(z)   at  z  =  0 

In  this  subsection  we  identify  and  briefly  discuss  the 
singularities  in  H(z) .  The  nature  and  origin  of  these 
singularities  are  the  same  as  in  Part  I,  and  the  discussion  given 
there  is  also  generally  applicable  to  the  present  case.  The  only 
difference  between  the  plane-strain  case  as  discussed  in  Part  I 
and   the  present   generalized  plane-strain   case   arises   because  in 
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the  former  case  the  integrand  for  H(z)  [see  eq  (54)  of  Part  I]  has 
E(q)  in  the  denominator,  whereas  in  the  present  case  it  has  the 
factor  (-Lq+0.5). 

The  singularities  in  H(z)  arise  from  (i)  the  homogeneous 
solution  and  (ii)  the  particular  solution.  These  are  discussed 
below: 

(i)  Singularities  in  the  Homogeneous  Solution 

The  form  of  the  singularity  is  given  by  eq  (104)  ,  where  krQ 
satisfies  the  constraint  given  by  eq  (97)  .  The  exponent  of  q  will 
have  a  value  between  0  and  -1.  This  term  will,  in  general,  have  a 
logarithmic     oscillatory     factor     unless  is      zero.  This 

singularity  will  be  present  only  if  a  Qr  [a  solution  of  eq  (76)  ] 
exists  which  satisfies  the  constraint  given  by  eq  (97) . 

(ii)  Singularities  in  the  Particular  Solution 

The  singularities  in  the  particular  solution  can  be 
classified  in  two  categories  according  to  the  nature  of  poles  in 
the  integrand  for  H(z)  as: 

(a)  Singularities  arising  from  simple  poles.  These  are  given  by 
the  first  term  in  eq  (81)  .  The  exponent  of  a  singularity  in  this 
category  can  have  a  value  between  0  and  -0.5.  These  terms  will 
also  have  an  oscillatory  log  factor  unless  Real   (Qr)   =  0. 

(b)  Singularities  arising  from  higher  order  poles.  These  are  given 
by  the  second  term  in  eq  (81)  .  There  is  at  least  one  second-order 
pole  at  q  =  -l/2  arising  from  the  factor  (-tq  +  0.5)  and  a  root  of 
D(q) ,  which  gives  rise  to  a  logarithmic  singularity  as  given  by  eq 
(86) .  This  singularity  has  no  oscillatory  factor. 
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If  D(q)  has  other  degenerate  roots,  H(z)  will  have 
singularities  of  the  form  given  by  eq  (87)  .  The  exponent  of  z  in 
this  class  of  singularity  may  have  a  value  between  0  and  -0.5.  The 
logarithmic  singularity  may,  in  principle,  have  a  power  of  up  to 
5.  There  will,  in  general,  be  a  logarithmic  oscillatory  factor 
unless  the  real  part  of  the  root  Qr  is  zero. 

The  exponents  of  all  the  singular  terms  are  independent  of 
the  loading  of  the  solid  and  depend  only  upon  the  elastic 
constants  and  the  geometry  of  the  solid.  This  is  because  the 
exponents  of  the  singularities  are  characteristics  of  the  matrix 
M,  which  does  not  depend  upon  the  type  of  the  loading.  The  loading 
only  affects  the  RHS  of  eq  (71) .  The  weights  of  the  singular  terms 
in  the  homogeneous  solution  are  also  independent  of  the  loading, 
as  is  apparent  from  eq  (75) .  However,  the  weights  of  the 
singularities  arising  in  the  particular  solution  will  depend  upon 
the  nature  of  the  loading. 

4.  Application  to  Z-5  Grain  Boundary  in  a  Cubic  Crystal 

In  this  section,  as  an  example  of  the  application  of  the 
formalism  developed  in  this  paper,  we  shall  apply  it  to  calculate 
the  stress  distribution  in  stainless  steel  containing  a  E-5  tilt 
grain  boundary.  We  consider  the  same  model  which  has  been 
described  in  section  6  of  Part  I.  We  consider  the  solid  subjected 
to  an  out-of -plane  load  such  that  e33  is  equal  to  a  constant  tj. 

First  we  consider  the  solid  to  be  infinite,  having  no  free 
surfaces.  The  stress  and  the  displacement  field  for  this  solid 
have  been  given  in  section  2.  The  expressions  for  various 
parameters  such  as  y  and  cr  have  been  given  in  section  6  of  Part  I. 

As  in  Part  I,  we  choose  the  units  such  that  c. .  =  1  for  solid 
A  in  the  UHP  and  tj  is  taken  to  be  unity.  The  calculated  values  of 
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ul(£)/  Tn(-)  '  T12  ^-^  anc*  T22^-^  can  ^e  eas^^y  obtained  from  eqs 
(34)-(39).  We  notice  from  these  equations  that  the  displacement 
field  is  a  linear  function  of  the  distance  from  the  origin,  and 
the  stresses  are  constants. 

The  calculated  values   of   the   components   of   the  vectors   T  . 

—ool 

and  as  defined  by  eqs    (36)-(39),    and  which  are  required  for 

eqs  (56)  -(59),  are  given  below  in  units  of  t)  [see  eq  (41)]: 

TA  ,   =   (-1.3,    0,    0)    and  TB  .   =   (-1.3,    0,  0) 

—  ool  —  ool 

and 

TA  0  =   (0,   -1.3,   0)   and  TB  0  =   (0,   -1.3,  0). 

—  03  2  —    00  2 

We  note   from  the  values  given  above  that  T  _    is  continuous 

-oo  2 

across  the  interface,  having  the  same  values  for  solids  A  and  B. 
This  is  required  by  the  boundary  condition  given  by  eq  (6)  .  In 
this    particular    case,     because    of    high    symmetry,  is  also 

continuous  and  equal  in  magnitude  to  T^ . 

Now  we  consider  the  effect  of  a  free  surface  in  this  model 
solid  using  the  formalism  presented  in  the  preceding  section.  The 
particular  solution  for  H(z)  is  given  by  eq  (70)  and  the 
homogeneous  solution  by  eq  (74)  .  In  order  to  see  whether  the 
homogeneous  solution  is  to  be  included  or  not,  we  have  to  solve  eq 
(76)  for  its  complex  roots  and  see  whether  it  has  any  roots  which 
satisfy  the  constraint  given  by  eq  (97) .  The  numerical  method  of 
calculation  of  the  roots  of  eq  (76)  has  been  described  in  the  next 
section.  The  two  roots  of  eq  (76)  which  have  a  positive  imaginary 
part  less  than  0.5  are 


Q  _  =  0.0272  +  0.4976  i  (105) 
^r0 


and 
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QrQ  =  0.0118  +  0.4994   i   .  (106) 

The   corresponding   singular   terms    in  H(z)    are   given   by  eqs 
(99)   and  (100) .  The  z-dependence  of  these  singular  terms  is 

H(Z)    *  z"0,9976cos(0.0272   In  z)  (107) 


and 


H(z)   «  z~0,9994cos(0.0118     In     z) .  (108) 


We  notice  from  eqs  (107) -(108)  that  these  singular  terms  in  the 
stress  also  have  a  logarithmic  oscillatory  factor  but  with  a  small 
magnitude. 

The  integrals  for  the  particular  solution  in  eqs  (61)  and 
(62)  are  calculated  numerically.  In  the  numerical  calculations  we 
have  taken  yw  to  be  10  and  the  limits  of  the  q — integral  from  -10 
to  +10.  For  testing  the  convergence  we  follow  the  procedure  given 
in  Appendix  I I -A. 

The  calculated  values  of  the  displacement  field  u1(x)  have 
been  shown  as  function  of  x.^  and  x2  in  a  3-D  graph  in  figure  II-2. 
We  notice  from  this  figure  that,  as  required  by  the  boundary 
condition  given  by  eq  (44),  the  u1(x)  is  continuous  across  the 
interface.  The  behavior  of  other  displacement  components  is 
similar  to  that  of  u^. 

The  calculated  values  of  the  stress  components  "t^fx)  and 
r22(x)  have  been  shown  as  functions  of  x.^  and  x2  in  3D  graphs  in 
Figs  II-3  and  II-4  respectively.  We  notice  from  these  figures 
that,  as  required  by  the  boundary  condition  given  by  eq  (45)  , 
these  stress  components  are  continuous  across  the  interface.  A 
similar  behavior  is  shown  by  the  stress  component  T22'  In 
addition,    the  stress  component  x12   is   zero  at  the  free  surface, 
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which  satisfies  the  boundary  condition  given  by  eqs  (46)  and  (47) . 
The  stress  component  t-^Cx)  is  shown  in  figure  II-5  as  a  function 
of  x^  and  x2 .  It  is  also  seen  to  be  zero  at  the  free  surface, 
which  satisfies  the  boundary  given  by  eqs  (46)  and  (47) .  It  is,  of 
course,  not  continuous  at  the  interface.  A  similar  behavior  is 
shown  by  the  stress  component  ^31» 

As  is  apparent  from  Figs  II-3  through  II-5,  the  stress 
components  are  singular  at  the  origin.  The  nature  of  the 
singularities  has  been  discussed  in  the  preceding  section.  The 
exponents  of  the  two  singular  terms  which  arise  from  the 
homogeneous  solution  are  determined  from  the  values  of  the  roots 
in  accordance  with  the  discussion  given  in  section  3D.  The 
behavior  of  these  singular  terms  is  given  by  eqs   (107)   and  (108) . 

In  order  to  identify  the  singularities  in  the  particular 
integral,  we  have  to  calculate  the  roots  of  eq  (76)  in  the  lower 
half  of  the  complex  q-plane.  As  remarked  in  the  previous  section, 
eq  (76)  has  one  root  at  q  =  -i/2.  The  corresponding  singularity  in 
the  H-function  will  therefore  be  a  log  singularity  as  given  by  eq 
(86) .  The  other  root  of  eq  (76)  with  negative  imaginary  parts  and 
which  satisfy  the  constraint  given  by  eq  (84)  is 


Q       =  0.0163   -  0.3342   i.  (109) 

pr  v  ' 

The  corresponding  singular  terms  in  H(z)   will  be  given  by  eq 
(81) .  The  nature  of  this  singularity  is 

H(z)   «  z~°*1658cos(0.0163  In  z) .  (110) 


We  find  that  this  singularity  also  has  a  logarithmic  oscillatory 
factor  with  a  small  coefficient.  As  remarked  in  section  3D,  the 
exponents   of   these   singularities   are   independent   of   the  loading 
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and  are  the  same  as  obtained  for  the  plane-strain  problem  in 
section  6  of  Part  I. 

5.  Application  to  Fiber-Reinforced  Laminated  Composites 

In  this  section  we  further  apply  the  formalism  developed  in 
this  paper  to  calculate  the  displacement  field  and  the  stress 
distribution  in  a  fiber-reinforced  laminated  composite  containing 
a  free  surface  and  subjected  to  an  out-of -plane  load.  We  shall 
consider  a  planar  interface  separating  two  different 
fiber-reinforced  layers,  terminating  at  a  stress-free  edge,  which 
is  taken  to  be  normal  to  the  interface.  As  mentioned  in  the 
introduction,  this  problem  is  important  for  studying  the  failure 
of  composite  materials. 

A  real  fiber-reinforced  laminated  composite  contains  several 
layers  and  interfaces  separating  them.  In  general,  the  interfaces 
are  not  planar.  However,  we  make  the  usual  simplifying  assumptions 
about  the  material  by  considering  it  as  a  bimaterial  composite 
with  a  single  planar  interface.  We  assume  that  both  the  layers  are 
homogeneous  and  anisotropic,  and  have  no  boundaries  except  for  the 
interface  separating  them  and  the  stress-free  edge.  In  spite  of 
various  approximations,  the  analytical  model  calculations  are  very 
useful  because  they  provide  a  good  understanding  of  the  underlying 
physical  principles  of  failure  of  composite  materials.  Moreover, 
they  provide  a  basis  for  an  efficient  numerical  approach  for 
calculation  of  the  stress  distribution  in  a  real  material  [14]. 

Thus  we  are  considering  a  model  in  which  the  two  solids,  A 
and  B,  in  figure  II-l  are  the  same  material  except  for  the 
orientation  of  the  fibers.  As  in  the  earlier  papers,  we  specify 
the  fiber  orientations  in  terms  of  0,  the  angle  of  rotation  about 
the  Y-axis  (x2~axis)  in  figure  II-l.  The  angle  9  is  measured  from 
the  Z-axis  as  shown  in  figure  II-6.   The  elastic  constants  of  the 
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two  solids,  A  and  B,  are  therefore  related  by  a  rotation 
transformation  given  by  the  orthonormal  matrix, 


S  = 


cos  0        0        sin  0  ] 
0  10 
-sin  0        0        cos  0 


(111) 


The  bimaterial  composite  in  figure  II-l  is  usually  specified 
in  terms  of  the  fiber  orientations  as  (8,0'),  which  means  that  the 
fibers  in  solids  A  and  B  are  at  angles  0  and  0',  respectively, 
from  the  Z  axis  as  shown  in  figure  II-6.  The  orientation  0=0 
refers  to  the  case  when  the  fibers  are  parallel  to  the  Z  axis. 

In  this  paper  we  shall  consider  composites  with  the 
orientations:  (0,90),  (30,-30)  and  (45,-45)  where  the  angles  are 
given  in  degrees.  There  is  a  strong  interest  in  the  engineering 
properties  of  composites  with  these  orientations.  We  shall 
consider  a  typical  high-modulus  graphite/epoxy  material.  The  set 
of  engineering  elastic  moduli  for  this  material  for  0=0,  as  given 
in  [10],  is 


■l  =  E2  =  1.447;   E3  =  13.79;   G12  =  =  =  0.586; 

"21  =  "31  =  "32  =  °-21 

where  E.   and  G. .   are  the  Young's  and  shear  moduli  respectively  in 
1        _  1  j 

units    of    10     kPa,    and  are   the    Poisson's    ratio.    The  above 

constants  are  based  upon  the  assumption  of  hexagonal  symmetry 
(transverse    isotropy)    and    are    the    same    as    used    in    [10].  The 

calculated  values  of  the  corresponding  elastic  constants,  which  we 

2 

shall  refer  to  as  Set  I,  are  given  below  (in  units  of  MN/m  ): 


Set  I 
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Cll  =  C22  =  15*25;   C12   =  °'329;   C13  =  C23  =  °'389; 
c33  =  139.5;  c44  =  c55  =  c66  =  0.586. 

In  addition  to  the  assumption  of  transverse  isotropy,  the 
elastic  constants  of  Set  I  are  based  upon  the  following 
assumption,  referred  to  as  the  Pipes-Pagano  approximation: 

G12  =  G31  and  V21  =  V32    *  (112i 

The  validity  of  the  assumption  in  eq   (112)   has  been  questioned  by 

Kriz    [13],    who   has   obtained   a   set   of   elastic   constants  without 

making  this  assumption.  This  set,  which  will  be  referred  to  as  Set 

2 

II,   is  given  below  in  units  of  MN/m    along  with  the  corresponding 

.    .         .  7 
set  of  engineering  elastic  moduli  in  units  of  10  kPa: 


Set  II 


Cll  =  C22  =  1'347;  C12  =  °'6747  C13  =  C23  =  °'618; 

c__  =  14.17;  c..  =  c,,  =  0.579;  oce.  =  0.336; 
33  44  55  66 

E1  =  E2  =  1.002;  E3  =  13.79;  G±2  =  0.336;  G23  =  G31  =  0.579; 

v21  =  0.491;  v21  =  v32  =  0.306. 

In  our  calculations,  we  have  included  both  the  sets  of 
elastic  constants  given  above.  The  calculations  are  carried  out 
exactly  along  the  lines  described  in  section  3  and  section  4.  The 
main  parameters  of  the  results  are  given  in  table  2 .  This  table 
gives  the  values  of  pa,  the  roots  of  the  Stroh  determinant,  as 
defined  by  eqs  (A. 19)  and  (A. 20)  of  Part  I,  for  a  =  1,  2  and  3  and 
for  solids  A  and  B.  These  roots  are  purely  imaginary  and  agree 
with  those  calculated  in  [11]. 
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TABLE  2 


Results  for  composite  materials  with  different  fiber  orientations 


Parameters 

Set  I 

Set  II 

(e,e')-> 

(0,90) 

(30,-30) 

(45,-45) 

(0,90) 

(30,-30) 

(45,-45) 

A 
Pi 

A 
P2 

A 
P3 

B 
Pi 

B 
P2 

B 
P3 

0.  8940L 
1.1185L 

1.  0l 

0.  6345L 

1.  0l 
4.7669L 

1. 0692L 
0 . 8626L 
2 . 53621 

Same 
above 

1. 0348L 
0.80931 
3 .44421 

Same 
above 

0.  96881 
1. 0322L 

1.  0L 

0 . 6820L 
1. 3126L 
4.75551 

1. 0915L 
0.9265L 
2 . 5306L 

Same 
above 

1. 1686L 
0.85461 
3.4362L 

Same 

as 
above 

TA  1 
col  2 
3 

0.  657 

0 

0 

4.617 
0 

6.873 

5.934 
0 

5.  301 

1.  067 

0 

0 

4  .883 
0 

6.  948 

6. 159 
0 

5.471 

TB  1 

ool  2 

3 

u .  /  uy 

0 

0 

a  cn 
4  .  bl  / 

0 

-6.873 

o .  y  j  4 
0 

-5. 301 

n  on 

0 
0 

4  .  883 
0 

-6. 948 

6. 159 
0 

-5.471 

TA,B  1 
oo  2  2 
3 

0 

0.614 
0 

0 

0.639 
0 

0 

0.613 
0 

0 

1.093 
0 

0 

1.078 
0 

0 

1.102 
0 

Roots  con1 

:ributing 

to  the  he 

Dmogeneous 

5  solutioi 

0. 00008 
0. 49900 

0 

0. 49992 
0 

-0.46563 
0 

-0. 49999 

I  Re 
Im 

II  Re 
Im 

Roots  con1 

0. 00175 
0.49978 

-0. 00011 
0.49993 

:ributina 

-0. 00351 
0.49522 

-0. 00097 
0.49905 

to  the  pj 

-0. 00588 
0.49812 

-0. 00019 
0.49981 

articular 

-0. 00181 
0. 49981 

0 

0.  49995 
solution 

-0. 00008 
0.49900 

0 

0.49992 
0 

-0.48102 
0 

-0.49999 

Re 
Im 

II  Re 
Im 

0 

-0.46662 

0. 00043 
-0.49928 

0 

-0.48839 

-0. 00069 
-0.49976 

0 

-0.47443 

-0. 00019 
-0.49981 

0 

-0. 38765 
0 

-0.49999 
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Further,  Table  2  gives  T  -  and  the  corresponding  values  of 
T  _  as  defined  in  eqs  (36)  -  (40)  .The  value  of  T  ,  is  needed  for 
calculating  the  RHS  of  eqs  (78)  and  (92)  .  We  immediately  notice 
that  both  the  stress  components  are  constants  as  shown  in  section 
2.  In  addition,  Tm2  has  the  same  value  in  UHP  and  LHP  as  required 
by  the  continuity  condition  given  by  eq  (6) . 

Table  2  also  gives  the  permissible  values  of  Q  in  the 
exponent  iQ-0.5  of  z  in  the  definition  of  H(z) .  These  values  of  Q 
are  the  roots  of  eq  (76) .  The  values  of  Q  which  contribute  to  the 
homogeneous  solution  given  by  eqs  (99)  and  (100)  are  those  which 
satisfy  the  constraint  given  by  eq  (97).  The  values  of  Q  which 
arise  in  the  particular  solution  are  those  which  satisfy  the 
constraint  given  by  eq  (84) . 

The  roots  arising  in  the  particular  solution  are  not 
explicitly  needed  for  the  calculation  of  the  stress  and  the 
displacement  field  since,  as  explained  in  the  text,  they  are 
included  in  the  integral  over  q  in  eq  (77)  .  However,  they 
represent  the  singularities  in  H(z)  and  the  stress  field  at  the 
origin  as  represented  by  eq  (81)   and  discussed  in  section  3D  (ii) . 

Table  2  shows  that  there  are  two  roots  which  contribute  to 
the  singularity  in  H(z)  as  described  in  section  3D(ii)  .  Out  of 
these,  the  exponent  of  the  singularity  arising  from  the  first  root 
agrees  with  that  obtained  by  Zwiers,  et  al.,  [10]  and  Ting  and 
Chou  [11]  for  the  same  composite.  This  root  has  no  real  part; 
therefore,  as  described  in  section  3  D(ii)  ,  it  will  have  no 
logarithmic  oscillatory  factor  in  agreement  with  [10]  and  [11]. 

The  second  root  is  very  close  to  -l/2.  As  described  in  the 
section  3,  it  can  be  shown  analytically  that  eq  (76)  has  a  root  at 
-l/2.  The  numerical  method  used  in  the  present  calculations  seems 
to  give  one  more  root  at  about  the  same  value  which  shows  that  the 
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root  at  -i/2   is  doubly  degenerate.    This  degeneracy  has  not  been 

reported  by  earlier  authors.   Since  this  root  of  eq   (76)    is  itself 

doubly  degenerate,  the  integrand  on  the  RHS  of  eq  (77)  will  have  a 

third  order  pole.  As  discussed  in  section  3  D(ii) ,  therefore,  this 

2 

implies  that  H(z)  has  a  [ln(z)]  singularity  at  z=0  and  not  just 
ln(z) .  The  weight  of  the  singularity  is  given  by  eq  (87) . 

The   roots  which   contribute  to  the  homogeneous   solution  are 

explicitly  introduced  in  the  total  solution  for  H(z) ,   as  given  by 

eqs    (99)    and    (100)  .    From  Table  2  we  see  that   in  addition  to  the 

singularity   in  the  particular   solution,    H(z)    will  have  two  more 

-1+e 

singular  terms   of   the   type   z  where   e    is   very   small — of  the 

-3  -5 
order  of  10       for  one  root  and  10       for  the  other  root.   The  two 

roots   are   apparently  degenerate,    which   indicates   a   higher  order 

pole  at  i/2.   However,   as  described  in  section  3,   the  root  at  l/2 

does  not  contribute  to  the  stress. 

The  complex  roots  of  eq  (76)  ,  as  quoted  in  Table  2  and  in 
section  4  for  the  cubic  case,  were  obtained  by  the  Muller's  method 
by  using  a  library  program  on  the  mainframe  computer  at  NIST. 
Equation  (76)  has  roots  at  ±l/2  which  are  known  analytically. 
These  were  given  as  known  roots  in  the  Muller's  method  and  were 
not  calculated  numerically. 

The  numerical  method,   however,   gave  two  apparently  new  roots 

close  to  this  value  as   shown  in  Table   2.    Since  these  roots  just 

satisfy  the  constraint  given  by  eq  (97) ,  they  were  included  in  the 

calculation.  We  could  not  ascertain  with  more  than  an  accuracy  of 
_  g 

10  that  these  are  independent  roots  (possibly  degenerate  with 
the  one  at  l/2)  or  the  same  root  at  l/2  which  is  being  picked  up 
by  the  numerical  method  with  a  round-off  error. 

In  order  to  verify  the  independence  of  these  roots,   we  first 

2 

divided  the  LHS  of  eq   (76)   by   (z  +0.25),   which  should  remove  the 
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known  roots  at  l/2  and  -i/2.   Then  we  applied  the  Muller's  method 

again   to   obtain   the   remaining   roots.    The  method   still   gave  the 

roots  as  given  in  the  table,   which  suggested,   within  an  error  of 
-5 

10  ,  that  these  are  indeed  independent  and  degenerate  roots. 
However,  the  uncertainty  is  not  serious  because,  as  shown  in 
section  3,  the  contribution  of  the  root  at  l/2  to  the  stress  is 
zero.  The  contribution  of  the  roots  in  Table  2,  which  are  close  to 
this  value,   is  also  quite  insignificant. 

We  calculated  the  integral  over  q  in  eq  (77)  numerically 
between  the  limits  -8  to  8.  The  convergence  of  this  integral  has 
been  discussed  in  Appendix  II-A.   The  value  of  y     was  taken  to  be 

10.    The    calculated    values    of    u^  T22    and    Tll    ^or  (°'9°) 

composite,  using  Set  I  elastic  constants,  have  been  shown  as 
functions  of  x.^  and  x2  as  3D  graphs  in  figure  II-7  through  figure 
11-10  respectively.  The  same  quantities  for  (30,-30)  and  (45,-45) 
composites  using  Set  I  elastic  constants  and  for  (0,90),  (30,-30) 
and  (45,-45)  composites  using  Set  II  elastic  constants  have  been 
shown  in  figure  11-11  through  11-3  0. 

We  have  given  the  above  results  for  composites  with  different 
orientations  because  of  a  strong  engineering  interest  in  these 
materials.  These  curves  should  be  useful  in  detailed  numerical 
studies  of  the  stress  distribution  in  real  composites. 

We  see  from  the  various  figures  that,  as  required  by  the 
boundary  conditions,  the  stress  components  t12  and  t22  and  the 
displacement  field  are  continuous  across  the  interface,  whereas, 
in  general,  is  not.  However,   in  case  of   (30/-30)   and  (45/-45) 

composites,  x11  is  also  continuous  across  the  interface.  This  is  a 
consequence  of  the  mirror  symmetry  across  the  interface. 
Apparently  for  the  same  reason,  u1  has  an  extremum  (maximum  or 
minimum)   at  the  interface  along  the  Y-  axis  in  these  materials. 
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Although  the  displacement  field  is  continuous  across  the 
interface,  in  case  of  the  (0/90)  composite  u.^  has  a  very  large 
gradient  at  the  interface  along  the  Y-axis.  This  is  a  consequence 
of  the  strong  elastic  anisotropy  of  the  solids  in  UHP  and  LHP  in 
the  composite.  This  can  be  understood  qualitatively  as  follows. 

The  stiffness  is  very  large  in  a  direction  parallel  to  the 
fibers  and  very  low  in  a  direction  perpendicular  to  the  fibers.  In 
case  of  the  (0/90)  composite,  the  fibers  in  the  UHP  and  the  LHP 
are  parallel  to  the  Z-  and  the  X-axes  respectively.  In  the  UHP  u1 
is  the  displacement  component  in  a  direction  perpendicular  to  the 
fibers,  which  is  against  relatively  weak  elastic  forces.  On  the 
other  hand,  in  the  LHP,  u.^  is  in  a  direction  parallel  to  the 
fibers,  which  is  against  relatively  strong  elastic  forces.  The  two 
are  therefore  very  different.  The  same  argument  explains  the  large 
difference  between  the  values   of  in  the  UHP   and  LHP   in  the 

(0/90)  composite. 

The    stress    components  and    t..    are    zero    at    the  free 

surface,  as  required  by  the  boundary  conditions  given  by  eqs  (46) 
and  (47) .  All  the  stresses  are  singular  at  the  origin,  as  shown  by 
the  peaks  in  the  figures.  These  peaks  appear  to  be  finite  in  the 
figures,  because  in  actual  calculations  we  calculate  the  stress 
close  to  the  origin  but  not  exactly  at  the  origin  for  obvious 
numerical  reasons. 

6.   Summary  of  Results  and  Discussion 

In  this  part  we  have  generalized  the  elastic  plane-strain 
Green's  function  method  given  in  Paper  I  and  used  the  method  of 
Part  I  to  calculate  the  stress  distribution  and  the  displacement 
field  in  a  bimaterial  composite  subjected  to  a  generalized 
plane-strain  loading.  The  stress  components  are  expressed  in  terms 
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of  a  function  H(z) ,  which  is  given  by  the  solution  of  a 
generalized  inhomogeneous  Hilbert  equation.  The  main  effort  is  in 
the  determination  of  the  function  H(z) .  The  final  result  contains 
two  parts:  a  particular  solution  which  is  expressed  in  the  form  of 
a  closed  integral  representation,  and  a  homogeneous  solution  which 
consists  of  a  finite  number  of  terms. 

The  integral  for  the  particular  solution  can  be  evaluated 
analytically  by  using  contour  integration.  This  gives  the  result 
in  the  form  of  a  series,  with  each  term  representing  the 
contribution  of  the  poles  in  the  integrand.  The  series  contain 
singular  as  well  as  nonsingular  terms  and  may  not  be  a  convergent 
series.  The  advantage  of  the  series  representation  is  that  it 
helps  in  precisely  identifying  the  singularities  in  the 
stress — their  weights  (the  coefficient  of  a  singular  term)  as  well 
as  exponents.  The  earlier  work  in  this  field,  in  general,  only 
attempted  to  identify  the  singular  terms  and  calculate  their 
exponents.  As  suggested  in  [10]  and  [11]/  the  weight  of  the 
singularity  has  to  be  obtained  by  a  numerical  solution  of  the 
boundary  value  problem. 

In  these  calculations  we  have  considered  only  one  free 
surface  in  the  solid  and  assumed  the  solid  to  extend  to  infinity 
in  other  directions.  The  effect  of  outer  boundaries  of  the  solid, 
which  in  real  cases  would  have  irregular  shapes,  is  very  difficult 
to  include  in  an  analytical  calculation.  This  can  be  easily  done 
by  using  a  numerical  technique  such  as  the  finite  element  method. 

However,  a  purely  numerical  technique  becomes  inefficient  and 
unreliable  in  a  region  where  the  solution  has  singularities  and 
discontinuities.  In  such  cases,  the  proper  procedure  would  be  to 
combine  the  analytical  and  the  numerical  techniques  by  using 
enriched  finite  element  method  [14]  or  taking  the  exact  analytical 
solution  near  the  origin  and  joining  it  by  numerical  solution  in 
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the  outer  regions.  The  formulation  given  in  the  present  paper 
should  be  particularly  suitable  for  this  purpose. 

In    this    paper    we    have    shown    that    the    stress    will  have 

—5 

singularities  of  the  form  z       (0  <  8  <  0.5)   and  ln(z)    (as  well  as 

—  1+5 

its  higher  powers)  arising  from  the  particular  solution  and  z 
arising  from  the  homogeneous  solution.  The  homogeneous  solution 
will  contribute  only  in  certain  cases  as  given  in  section  3  D  (i) . 
In  addition,  depending  upon  the  material  constants,  the 
singularity  may  have  a  logarithmic  oscillatory  factor  which  is 
similar  to  that  associated  with  interfacial  cracks  (see,  for 
example,  Paper  II) . 

We  did  not  find  the  logarithmic  oscillatory  factor  in  case  of 
the  laminated  composites  which  is  in  agreement  with  the  earlier 
work  on  these  materials.  However,  in  case  of  the  cubic  solids 
described  in  the  preceding  section,  we  did  find  this  factor, 
although  with  a  small  magnitude.  Our  result  on  the  nature  of  the 
singularities  of  the  type  z  5  and  ln(z)  in  laminated  composites 
agrees  with  that  obtained  in  some  earlier  papers  [10,11].  However, 
the  singularities  with  exponents  between  -0.5  and  -1.0  have  not 
been  identified  in  the  earlier  papers. 

The  existence  of  the  homogeneous  solution  in  our  result  is 
due  to  a  boundary  condition  which  we  impose  on  the  solution  of  the 
Hilbert  eguation.  This  condition  requires  the  displacements  to  be 
single-valued.  It  is  obviously  a  physical  requirement  which  comes 
naturally  in  our  method  of  solution  because  otherwise  the  Hilbert 
equation  will  not  have  a  unique  solution.  This  seems  to  have  been 
missed  in  the  earlier  published  work.  This  may  be  the  reason  why 
the  singularities  arising  in  the  homogeneous  solution  with 
exponents  between  -0.5  and  -1.0  have  not  been  identified  so  far. 

It   can   be   easily   seen   analytically  that   if   the   exponent   of  the 
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singularity  in  the  homogeneous  solution  is  exactly  -1,  it  will  not 
contribute  to  the  stress.  In  the  case  of  the  laminated  composites 
which  was  considered  in  section  5,  the  only  terms  which  we  found 
in  the  homogeneous  solution  have  exponents  very  close  to  -1 
(imaginary  part  of  Q  close  to  0.5).  These  terms  do  not  make  a 
significant  contribution  to  the  stress  and  therefore  do  not 
invalidate  the  results  published  in  the  earlier  papers  on  such 
laminated  composites.  Even  in  the  cubic  case  which  was  discussed 
in  section  4,  the  exponents  of  these  terms  are  quite  close  to  -1. 

Regarding  the  singularities  in  the  particular  solution,  we 
seem  to  obtain  a  degenerate  root  at  Q  =  -i/2  as  shown  in  Table  2. 

As  remarked  in  the  preceding  section,   this  implies  a  singularity 

2  ...  .  .  . 

of   the   type    [In   z]    .    The  possibility   of   such   singularities  has 

been  acknowledged  in  [10]  but  not  actually  included  in  the  earlier 

calculations.     In    our    method,    we    do    not    need    to    isolate  the 

singularities      in      the      particular      solution     which      are  all 

automatically  included  in  the  integral  representation. 

One  of  the  advantages  of  our  method  is  that  we  do  not  have  to 
calculate  the  exact  roots  of  eq  (76)  except  for  the  homogeneous 
solution  when  it  exists.  In  these  cases,  we  have  only  to  search 
for  the  roots  in  a  finite  range  as  given  by  eq  (97)  .  Moreover, 
even  when  we  have  to  calculate  the  roots,  we  need  to  solve  only  a 
6x6  determinantal  equation.  In  the  earlier  methods,  [10,11]  one 
has  to  solve  a  12  x  12  determinantal  equation.  The  main 
disadvantage  of  this  method  is  that  it  cannot  account  for  the 
other  surfaces  and  interfaces  which  exist  in  a  real  solid.  In  real 
cases,  a  combination  of  a  numerical  and  the  present  technique,  as 
suggested  earlier,  should  provide  an  efficient  method  for 
analyzing  the  stress  distribution. 
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Appendix  I I -A 


Transform  of  a  Constant  and  Some  Numerical  Considerations 

In  this  appendix,  we  describe  the  transform  of  a  constant 
which  is  required  for  calculating  the  RHS  of  eq  (63)  .  In  this 
context  we  shall  also  discuss  some  aspects  of  the  numerical 
convergence  which  are  relevant  to  the  calculations  in  this  paper. 

Let  C  denote  a  constant,  which  is  expressed  in  terms  of  its 
transform  c(q)  as 


00 

C  = 

—or 


c(q)yLq  °*5dq  .  (A.l) 


Following  the  method  given  in  Appendix  I-B  of  Part  I,  we  obtain 
the  following  inverse  relation  for  c(q) 


C  y"1-^0-5      ,  (A. 2) 


0J 

which  is  an  equation  of  the  same  form  as  eq  (64)   or  eq  (65) . 

We  evaluate  the  integral  in  eq   (A.l)    in  the  limits  0  and  y 


00 


and  take  the   limit  as  y    goes  to   infinity.    The   integral   is  then 

00 

elementary,  and  we  obtain  the  following  expression  for  c(q) 


*   x  CI  ,     ,-Lq+0.5  -i 

C«3)   =  "  JWZ   (q+0.5L)  *  (A*3) 

In  order  to  verify  that  this  gives  the  correct  transform  of  the 
constant  C,  we  substitute  for  c(q)  from  eq  (A. 3)  in  eq  (A.l). 
Obviously   the    inverse   relation    is   valid    only    if    the  following 
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integral  is  independent  of  y  in  the  limit  y  -»  oo  and  is  equal  to 
unity;  that  is, 


00 


27TL 

-00' 


(q+0.5L)  woo 


(YjY)    Lq+°'5dq  =  1   .  (A. 4) 


The  integral  on  the  LHS  of  eq  (A.  4)  can  be  easily  evaluated 
by  taking  a  semicircular  contour  in  the  lower  half  plane  as  shown 
in  figure  1-4  of  Part  I.  In  the  limit  of  large  Ym  the  contribution 
of  the  semicircle  will  vanish,  and  the  integral  can  be  seen  to  be 
unity  as  required  by  eq  (A. 4)   for  all  values  of  y. 

Equation  (A. 4)  can  be  used  for  testing  the  convergence  of  the 
numerical  integration  over  q  in  eqs  (61)  and  (62)  .  In  order  to 
carry  out  the  numerical  integration  for  H(z)  ,  we  have  to  choose 
finite  but  large  values  of  y^  and  the  limits  of  the  integration 
for  q.  We  also  evaluate  numerically  the  integral  in  eq  (A. 4)  with 
the  same  values  and  check  whether  the  integral  is  acceptably  close 
to  unity. 

One  problem  which  arises  in  this  context  is  that  if  the 
limits  of  integration  for  q  are  taken  to  be  too  large,  the 
integrand  in  eq  (61)  may  become  numerically  unstable.  Some  matrix 
elements  of  the  matrix  M(q)  contain  factors  like  exp(±7rq),  which 
will  become  too  large  near  the  upper  or  the  lower  limit  of  q  in 
the  integration.  The  matrix  will  then  become  ill-conditioned  and 
its  calculated  inverse  may  contain  large  numerical  errors.  On  the 
other  hand  if  the  limits  of  q  are  taken  to  be  too  small,  the 
integral  may  not  converge  in  the  sense  required  by  eq  (A. 4). 


The  stress  T1(x)  as  calculated  from  eqs  (58)  and  (59)  will  be 
zero  at  x1  =  0  only  if  eq  (A.  4)  is  satisfied.  It  is  therefore 
necessary  to  choose  optimum  values  for  the  limits  of  q  so  that  eq 
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(A. 4)  is  satisfied  without  making  the  matrix  M(q)  ill-conditioned. 
To  make  the  matters  worse,  the  convergence  also  depends  upon 
values  of  y  in  eq  (A. 4).  This  would  suggest  that  different  limits 
for  the  q-integral  have  to  be  chosen  for  different  values  of  y. 
This  is  not  advisable  because  it  might  introduce  artificial 
discontinuities  in  the  result. 

One  way  out  of  this  difficulty,  which  we  have  used  in  our 
calculations,  is  the  following.  We  choose  reasonable  values  for  y^ 
and  the  limits  of  q  and  evaluate  the  integral  over  q  for  H(z)  in 
eqs  (61)  and  (62)  using  these  constants.  Then  we  calculate  the 
integral  in  eq  (A. 4)  by  using  the  same  values  of  y,  y  ,  limits  of 
q  and  the  integration  interval.  The  calculated  value  of  the 
integral  in  eq  (A. 4)  is  taken  to  be  the  effective  value  of  unity. 
We  then  divide  H(z)  by  this  effective  value  of  unity.  The  errors 
due  to  lack  of  convergence  are  then  largely  cancelled  out  from  the 
result.  This  is  indicated  by  the  fact  that  T^x)  as  calculated 
from  eqs  (58)  and  (59)  will  become  zero  at  X-  =  0  in  accordance 
with  the  free-surface  boundary  condition. 
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fB(x)  ////// 


Solid  A  (UHP) 


////// 


fw  \\\\\\° 


www 


Solid  A  (LHP) 


Figure  II-l:     An  infinite  bimaterial  composite  containing  a  plane  interface 
and  the  coordinate  system  used  in  these  calculations.     The  force  functions 
f  (x)  and  f  (x)  are  applied  just  outside  the  UHP  and  LHP  respectively. 
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Figure  II -6:     The  rotation  angle  of  the  fibers  with  respect  to  the  Z  or  the 
X  axis.     The  rotation  is  about  the  Y  axis  which  is  normal  to  the  plane  of 
the  paper.     At  0  orientation,  the  fibers  are  parallel  to  the  Z  axis.  At 
orientation  9,   the  fibers  are  parallel  to  OZ' . 
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It  is  found  that,  as  predicted  by  earlier  authors,  the  stress  is  singular  at  the  intersection  of  the  free  surface 
sad  the  interface  in  both  plane  strain  and  generalized  plane-strain  case.     These  singularities  in  the  stress  field 
associated  with  the  presence  of  the  fret  surface  are  identified  and  discussed.     The  singularies  in  the  stress  field 
are  found  to  be  of  the  type  r~& ,  ln(r)  as  well  as  containing  higher  powers  of  ln(r),  where  6  is  between  0  and  1  and 
r  is  the  radial  distance  from  the  intersection  of  the  free  surface  and  the  interface.     In  addition,  it  is  found 
that,  in  general,  the  stress  field  will  also  contain  an  oscillatory  factor  of  the  type  exp  [  (.gln(  r)  ] ,  where  g  is  a 
constant  that  depends  upon  the  material  parameters  cf  the  two  solids.     However,  for  the  case  of  the  generalized 
plane  strain  in  fiber-reinforced  composites,   the  factor    g    is  almost    0    so  that  the  oscillatory  behavior  is  absent. 

Part  I  also  gives  a  simple  and  convenient  method  for  solving  an  inhomogeneous  generalized  vector  Hilbert  problem, 
which  contains  a  nonsingular  kernel  in  addition  to  the  usual  singular  Hilbert  kernel.     The  solution  of  this  Hilbert 
problem  is  required  for  obtaining  the  Green's  function  in  the  present  case,  as  well  as  many  other  cases  dealing  with 
the  stress  analysis  of  composite  materials  containing  surfaces,  interfaces,  and  cracks. 

The  solution  of  the  Hilbert  problem  is  obtained  by  using  a  complex  transform  of  the  type  y ^  where  y  and  q  are 

variables  on  the  real  axis.  It  is  shown  that  this  function  is  an  eigenf unction  of  the  Hilbert  operator  and  is  ortho- 
gonal for  y,  varying  in  the  range  between  0  and  *  oo    and  q  between  -<»  to  +  <». 
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